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ABSTRACT 


\Ql9°^ 

Digital (numerical) filtering techniques have became significant 
methods for data processing. This report presents the necessary back- 
ground theory for a large class of digital filters. Particular filters 
discussed are those for smoothing, smoothing and first derivative, 
smoothing and second derivative, and smoothing and integration. Error 
bounds are given for the first three types. 




PREFACE 


Under Contract No. NAS8-5164 with NASA, Auburn Research Foundation 
undertook a program to study numerical smoothing and differentiation 
methods and numerical aspects of finite difference methods. The work 
was performed by the Mathematics Department, Auburn University* The 
Senior Investigators were Dr. Nathaniel Macon, 1 October 1962 to 1 June 
1963, and Dr. L. P. Burton, 1 June 1963 to 31 July 1964- . Investigators 
were Edward B. Anders, 1 October 1962 to 1 June 1963, Paul W. Spikes, 

1 October 1962 to 1 September 1963, Alfred D. Lasaine, 1 October 1962- 
to 1 June 1964, James T. Taylo, 1 June 1963 to 20 July 1964, and James 
J. Johnson, 1 September 1963 to 1 June 1964. 

The main effort was devoted to linear digital (numerical) filters 
for performing the smoothing, differentiation, and integration of 
discrete data and to error analysis for these filters. The primary 
interest was filtering techniques for one variable, but an extension 
to n variables was developed by Anders [l4]. 

This report incorporates the necessary background theory starting 
with the classical Fourier theory and going into generalized functions, 
pertinent results obtained by other workers, and the results and con- 
clusions obtained under this contract. Numerous references to other 
publications are given. 

A reader interested only in the application of the techniques 
discussed here may start reading with Chapter IV. 

E.B.A. 

* J. J. J • 

With Ronald J. Graham and David G. Aichele A.D.L. 

as NASA Technical Representatives. P.W.S. 

J.T.T. 
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CHAPTER I 


CLASSICAL FOURIER ANALYSIS 


In this chapter we shall list, without proof, some of the 
results from Classical Fourier Analysis. Proofs of these restats are 
readily available in many good texts on the subject [1], [2]. 


Fourier integral theorems 

Given a function f(x) of a real variable x, if f(x) is absolutely 
int egr able , 1 . e . , if 


oo 


/ |f(x) |dx < 00 , 


then the integral 


exists and 


oo 

F(y) = f f(x) 
J -oo 

00 

f ( x) = / F(y) 

-CD 


-2jtixy , 
e dx 


2nixy. 
e dy. 


( 1 . 1 ) 


( 1 . 2 ) 


(1.3) 


The condition (l.l) is sufficient but not necessary for (1.2) 
and ( 1. 3) . We assume that all functions with which we shall be 
concerned in this chapter satisfy conditions necessary for the 
existence of (1.2) and (1.3). 

The usual notation 

f(x) < > F(y) 

will be used to indicate that the functions of f(x) and F(y) are related 
by (1.2) and (l. 3). We call F(y) the Fourier transform of f(x) , and 
f(x) the inverse Fourier transform of F(y) . This being the only type 
of transform we shall use, it is occasionally convenient to drop the 
word "Fourier" and speak of the "transform of f(x)" and of the "inverse 
transform of F(y)." 



(1.4a) 


F(y) is in general complex: 

F(y) = R(y) + il(y) 
or 

F(y) = A(y)e i0(y) 

2 

A(y) is called the Fourier spectrum of f(x), A (y) is its ener 
spectrum , and e(y) its phase angle . 

The function f(x) may in general be complex: 

f ( x) ■ f R (x) + ifj(x) 

and „ 

00 

F(y) = / [f_(x) + if T (x) ] [cos 2:txy - i sin 2*xy]<3x 
R 1 


(1.4b) 


UU 

if [f (x) cos 2nxy + f_(x) sin 2itxy]dx 


- if [f (x) sin 2-rtxy - fj(x) cos 27txy]dx. 
-00 


Hence 


(1.5) 


r°° 

R(y) = / [f (x) cos 2itxy + f-^x) sin 2«xy]dx ( 1.6) 


r oo 

[f (x) sin 2rtxy - f-^x) cos 27txy]dx. (1.7) 


Using (l.4a) we find, in a similar manner, that 


oo 

j( x) = f [R(y) cos 2nxy - l(y) sin 27txy]dy. (l.8) 


r 

f (x) * I [R(y) sin 2«xy + l(y) cos 2«xy]dy. (l-9) 
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I 


Suppose that f(x) is a real function of x. Then 


oo 

R(y) = f f(x) cos 2jtxydx 


and 


-oo 
oo 


i(y) - f f ( x) sin 2itxyclx. 


-oo 


It is easily seen that 


R(-y) = R(y) 


and 


i(-y) = -i(y), 

i.e. , R(y) is even and l(y) is odd. 

Hence 

F(-y) = F*(y) 

where the asterisk indicates the complex conjugate. 
Conversely, if (1.12) holds, then 

oo 

f T ( x ) = / fH(y) sin 2itxy + l(y) cos 2«xy]dy 

-oo 

■ 0 

because the integral is odd, and therefore f(x) is real. 

Suppose that f(x) is an imaginary function of x. Then 

oo 


R(y) f( x) sin 2nxydx 


and 


-oo 

oo 


l(y) = / f(x) cos 2nxydx. 

-oo 

Hence R(y) is odd and l(y) is even, and 

F( -y) = -F*(y). 

Conversely, if (l.lU) holds, then f(x) is imaginary. 


( 1 . 10 ) 


( 1 . 11 ) 


( 1 . 12 ) 


(1.13) 


(1.1V) 
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Even and odd functions 


If F(y) is even, then 
r .oo 

f(x) = 2/ F(y) cos 2«xydy. 
o 

And if F(y) is odd, then 
oo 

f(x) = sin 2,txydy ' 


(1.15) 


(1.16) 


Similar expressions are found for F(y), f(x) being even or odd. 

We will now list some of the more important theorems of class! 
cal Fourier analysis. 


I. Linearity 

If F(v) and G(y) are the transforms of f(x) and g(x) 
respectively, and if a,b are arbitrary constants, then 

a F(y) + b G(y) < > a f(x) + b g(x). (1.17) 


II. Symmetry 

If F(y) is the transform of f(x), then 

F(x) < » f( -y) • (1-18) 

If a is a non-zero real constant and F(y) is the transform 
of f ( x) , then f(ax) < > jjj- F(£j. (1.19) 


IV. "x" domain shifting 

If x Q 
f(x) , then 


jf Xq is a real constant and F(y) is the transform of 


f(x-xo) < 5* F(y) e‘ 2niX o y 


( 1 . 20 ) 


V. 


"y" domain shifting 

If y 0 is a real constant and F(y) is the transform of 
f(x) , then 


e 2mxy 0 f ( x ) < > F(y-y c ) 


( 1 . 21 ) 
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Note: 


Using (1.19) and (l.2l) we have 


2jtixy 0 


f(ax) < (1-22) 


VI. "x" domain differentiation 


If the transform of d n f(x ) exists and if F(y) is the 


dx 


n 


transform of f(x), then 


d n f(x) « > (23tiy) n F(y) 


dx 


n 


(1.23) 


VII. "y" domain different iat ion 

If the inverse transform of d I V(y ) exists and if F(y) 


<Jy 


n 


is the transform of f(x), then. 

(-2*ix) n f(x) <€ > d n F(y) 


(1.24) 


n 


cly 

VIII. Conjugate functions 

If F(y) is the transform of f(x) and the asterisk 
indicates the conjugate, then 

f*(x) < > F*( -y) . 

IX. "x" domain convolution 


(1.25) 


If F(y) is the transform of f(x) and G(y) is the transform 
of g(x), then 


OO 

f(x)*g(x) = f f(z)g(x -z)dz < > 


F(y)G(y) 


(1.26) 


•00 


We note here that the existence of the transforms of f(x) 
and g(x) is not sufficient to prove this theorem. The proof 
usually goes as follows: 


Let 


/ CD 00 

p -2«ixy T f( 
-00 ^00 


CD 

f( z)g(x-z)dz 


dx. 
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Assuming that the order of integration can be changed, we 


have 


oo 

H(y) = F f(z) 
^-00 

Using (1.20) we have 


i 


00 

P -2jtixy t \ , 

e g^x-z;ax 

oo 


az« 


00 

H(y) = F f(z)e” 2rfiyZ G(y)dz 

vim 


^-00 
*.F(y)G(y). 

A sufficient condition which would allow the inter- 
change in the order of integration is that f(x) and g(x) 
be square- integrable, i.e. 

oo 


J |f(x)| ! 


dx < oo 


and 


-oo 
oo 

/' g( 


x) I dx < oo* 


(1.2T) 


Similarly, we have the following theorem which can 
be proved from ( 1 . 26 ) by using (l.l8). 

X. "y" domain convolution 

If F(y) is the transform of f(x) and G(y) is the 

transform of g(x), then 

oo 

f(x)g(x) < ^ F F(z)G(y-z)dz = F(y)iK}(y). (1.28) 

J- oo 

XI. Parseval* s formula 

If F(y) and G(y) are the transforms of f(x) and g(x) 
respectively, then 

OO ^00 

(1.29) 


w w 

F f(x)g(x)dx * / F( -y)G(y) dy. 

^ -OO -00 
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CHAPTER II 


GENERALIZED FUNCTIONS AND THEIR FOURIER TRANSFORMS 


2.1 Generalized functions 

The Dirac delta functions (t) is often "defined" by one of the 
following statements: 


(A) 


If f(t) is a continuous function at t=t Q , 
the property that 

GO 



*(t-t Q )dt = f(t 0 ) ; 


then S (t) has 


( 2 . 1 ) 


(B) S( t) = o if t^o, and 


oo 


J $(t)dt = 1 ; 


( 2 . 2 ) 


-oo 

(C) S(t) * lim f n (t) where |f n ( t )"l is a sequence of 
n-*oo 1 J 

functions satisfying the conditions: 

r 00 

/ f n (t)dt = 1 and lim fn(t) = o, t^o. 

-oo 


(2.3) 


n->oo 


These definitions are meaningless if we attempt to think, of &(t) 
as a function in the ordinary sense. By introducing the delta function 
as a new concept, a generalized function, the difficulties are 
resolved. Then (.2.1) can be given precise meaning, but (2.2) and 
(2.3) do not uniquely describe S(t). Such an approach to the pro- 
blem is usually relegated to an appendix of a book, and is usually 
very sketchy. 

We are concerned here with extending the Fourier integral theory. 

A relatively easy and short approach, which is similar to the extension 
of the rational numbers to the real numbers, is given by Lighthill [ 3 ]. 
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Alternately, we could delve into functional analysis to extend the 
theory. For this study it appears that the first approach will be 
more meaningful, and so we have chosen to develop it. However we 
shall restrict ourselves to that portion of the development which 
suit s our purpose here, and for a detailed development the reader 
is referred to Lighthill. 

Advantages in clarity are gained by using some of the term- 
inology and notation of functional analysis. 

Let R be a set of scalars (the real or complex numbers) and let 
K be a set with an operation denoted by (+) and called "addition" 
defined on it. 


Definition 1 . K is called a linear space over R if 

1) K is an abelian group with respect to (+) 

2) ak e K for all a e R and all k e K, and 

(a) a(bk) = (ab)k for all a,b e R, 

(b) 1 • k = k 

3) If a,b e R and k , k 2 € K, then 

(a) a(k x +k 2 ) = ak^+a^ 

(b) (a+b)k^ = ak^+bk^. 

Let 

X = |x| x is a real numbei^ 

and 

a = |h| h is a function of x, x a real number 

Under ordinary addition and scalar multiplication, H is a linear 
space over X. 

We denote the Fourier transform of f(x) by Z(f), 


Z(f) = 



-2flixy , 

dx. 


We denote the inverse transform by Z (g) = f. Then Z is a mapping 

-1 -1 
of a subset Z of H into H, and Z is a mapping of a subset Z 

X X 
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-1 

of H into H. However Z =j= Z . that is, the domains of definition of Z 
and Z in H are not identical. In addition, there exist convergent 
sequences of functions of Z x , each element of which has a transform, 
but whose limit has no transform. 

Z and Z “ 1 

x x are linear spaces. This suggests the possibility of 

constructing a linear space in which these may be imbedded and such 
that the mappings Z and Z _1 can be suitably extended to remove the 
above difficulties. 

Let S denote the set of all functions f(x) which are everywhere 
differentiable any number of times and such that f(x) and all of 
its derivatives are 0( jx|~ N ) as |x| — » oo for all integers N. As 
a reminder, the 0 notation, f(x) = 0(g(x)) as x — > a, means that 
there exists a positive constant A such that 

1 f(x) | < A| g(x) | 

as x — > a. g 

Example: e" x is contained in S. 

It is clear that S, under ordinary addition and scalar multi- 
plication of functions, is a linear space. S is called a test or 
fundamental space and elements of S are called test functions . 

If g(x) is a function such that g(x)f(x) € S for all f(x) € S, then 
g(x) is called a multiplier on S. 

Let M denote the set of all functions m(x) which are everywhere 

differentiable any number of times and such that m(x) and an of its 

N 

derivatives are 0( |x| °) as |x| — > oo for some integer N q . 

Lemma : If m(x) e M and f(x) e S, then m(x)f(x) e S. 

It suffices to show each term on the right in the above equation 
is in S. There exist numbers A^ > 0, > 0 and an integer 

such that 

/.x N 

|m'' , '(x)| < A 1 |x| for all x such that |x| > K^. 
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If N is any integer, then there exist numbers A > 0, K > 0 such 


|f Vp " d; (x)| < Ag |x| 2 for all x such that |x| > Kg. 

Then for all x such that |x| > max IX^Kg] , 


|m (j) (x)f (p ' d) (x) | <A x A 2 |x| 1 2 . 

But N = N^-Ng is arbitrary since Ng is arbitrary. Hence m^(x) 

f(p-o)( x ) e s. This shows that the elements of M are multipliers on 
S. Obviously if m(x) € M, then m' (x) and ra'(x)f(x) e S. 

Example: Any polynomial is contained in M. 

Theorem 1 . If f(x) e S, then 

(a) f'(x) e S 

(b) Z(f) = g(y) e S 

(c) Z -1 (f) = h(y) € S 

* (d) f(-x) € S 

(e) f*(x) e S 

(f) f(ax+b) e S, a,b constants and a ^ 0. 

Proof: The first part is obvious. To prove (b), let 


g(y) = f f(x) e” 2rtixy dx 


and differentiate p times. Then integrate by parts N times « 


g (p, (y> 


00 N 

^”N / ~N {(-2«ix) P f(x)} 

(2«iy) N J .oo dx N 


< r f x p f ( x i 

|y| ^00 , N l x f( } 
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Therefore 


g^(y) = 0( |y|” N ). 

Part (c) is proved in a similar manner. For part (d), from (b) 
we have g(y) = Z(f) is a test function. Applying Z 1 we have 


oo 

f(x) = f g(y) 
J-CO 


2*ixy , 
e dy. 


Hence 


oo 


f(-x) = f g(y)e" 2rtixy dy, 
-oo 


i.e., f(-x) = Z(g), hence f(-x) € S. 

Part (e) follows by virtue of the fact that |f^ m ^(x)| = |f*( m ^(x) | . 
For part (d), let f(x) < > g(y) and use (1.19) and (l.20) to write 


f(ax+b) <- 


a 


«(?). 


(a. 5) 


Then it suffices to show that the right side of (2.5) is a test 
function. Clearly it is everywhere differentiable any number of 
times. Since 

2«i— ^ /jr\ 

it is obviously 0(|y|' N ) as |y| — > oo for all integers N. Hence 
f(ax+b) e S. 

Definition 2. A sequence of test functions is called 

regular if 

oo 

lijn ^ J f n (x)F(x)dx < oo for all F(x) e 8. (2.6) 

\i] 

An example of such sequence is ^e n J • 

Definition 3 . Let C denote the class of all regular sequences. 

Then {f^U)} e C is said to be equivalent to {g^xjj. e C if 
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oo oo 

/ f (x)F(x)dx = lim / g (x)F(x)dx for all F(x) € S. (2.7) 
J ~ n 00 n 


It is easy to see that this is an equivalence relation. Hence 
C is partitioned into disjoint subclasses, and we let S denote the 
set of all subclasses of C determined by this equivalence relation. 

Definition 4 . An element of S is called a generalized function . 

Hence a generalized function is the class of all regular sequences 
equivalent to a given regular sequence. Since the limit (2.6) is 
the same for all sequences of a given element s of S, any sequence 
of s can serve as a representative of that class. 

Definition 5 . If s(x) € S and F(x) e S we define 


co 00 

I s(x)F(x)dx = U m / s (x)F(x)dx 


( 2 . 8 ) 


where { s n (x) } is any representative of s(x). 

In general, the left side of (2.8) has no meaning as an integral 
in the ordinary sense. It is just a conventional notation we adopt 
to denote the limit on the right side of (2.8). 


Example: The sequence Ip 11 Jrepresents a generalized function l(x) 


such that 


p 00 p 00 

/ l(x)F(x)dx = / F(x)dx. 
'■too -00 


(2.9) 


Definition 6 . Let f(x) and h(x) be generalized functions with repre- 
sentative sequences and {^ n (x)} respectively. Also let 

m(x) e M. Then 

(a) the sum f(.x) + h(x) is defined as the generalized function 
b(x) with a representative {f n ( x) * h n W} ; 

(b) the derivative f'(x) is defined as the generalized function 
with a representative te(*>} ; 
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and 


(c) f(ax+b) is defined as the generalized function with a 
representative {f^(ax+b)| ; 

(d) The product m(x)f(x) is defined as the generalized function 
with a representative {ra(x)f n (x)j’ ; 

(e) The Fourier transform g(y) of f(x) is defined as the 

generalized function with a representative ^g (y) ^ y 

where g (y) is the Fourier transform of f (x). 
n n 

To show that these definitions are consistent we must show that 

(a) each sequence named is a sequence of test functions, 

(b) each sequence named is regular, 

(c) equivalent regular sequences defining f(x) ; and h(x) 
lead to equivalent regular sequences defining the new 
generalized functions. 

Part (a) follows in each case from previous remarks and Theorem 1. 
Parts (b) and (c) can be deduced in each case from the following 
equations: 

r°° r°° 

11 h^koj [f „ (x)+h n (x)lF(x)d * * f n ( x,F(x,dJC 

-00 -oo 

00 

+ &ooJ h n ( x ) F ( x ) dx » (2.10) 

-oo 

00 oo 

2 ) &oo J f A (x)F(x)dx = -li^ jTf^( x )F’( x )dx, (2.11) 

CD OO 

3) *^A (ax+b)F(x)dx ^oo A (x)F(2 ? )dx (2 - i2) 

-oo -oo 

00 

*0 IjyjJoo I [m(x)f (x)]F(x)dx = 

-00 

r°° 

n^J f n ( x )[m(x)F(x)]dx (2.13) 
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and, from Parseval's formula. 


n <X> oo 

^ nr^oo /h (y)G(y)dy = lijn f f (x)F(-x)dx, 
J -oo n n ^°° X.qo“ 


(2.14) 


where G(y) is the Fourier transform of F(x). 

It follows from the equations (2.10) — (2.l4) and definition 5 

that 


and 


oo oo 


/ f’(x)F(x)dx 

= - f(x)F'(x)dx, 

(2.15) 

-00 

—00 


OO 

n 

r°° 


/ f(ax+b)F(x)dx 
-00 

=Rj f(x)r(^)ax ( 

‘ 1 -OO 

(2.16) 

00 

n 

OO 


/ [m(x)f(x) ]F(x)dx 

= / f(x)[m(x)F(x)]dx, 

(2.17) 

-CX) 

-00 


00 

n 

00 

p 


/ My)G(y)dy 

= / f(x)F(-x)dx 

(2.18) 

-00 

-00 



If a is a scalar, then the function h(x) = a is contained 

in M and hence if f(x) e &, af(x) e S. With the addition as defined 

above in S, we have S is a linear space. 

Let I denote the set of all ordinary functions h(x) such 
2 -N 

that (l+x ) h(x) is absolutely integrable on (-oo, oo) for some 
integer N. 


Theorem 2 . I is a linear space and can be imbedded in S, i.e., ICS. 

The theorem asserts that if h(x) e I then there exists a regular 
sequence fn (x)} such that 


r°° r°° 

lijf^ J h (x)F(x)dx = / h(x)F(x)dx, for all F(x) € S. (2.19) 
-oo -00 


For a proof of this theorem we refer the reader to [3], pp. 22-23. 
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Having shown the existence of such a sequence, the mapping 

** C { h n} > 

where c|h n |denotes the class of. all sequences equivalent to 
^h n j. , imbeds I in S. 

Note that the integral on the right side of (2.19) is 


oo oo 

J h(x)F(x)dx = J {(l + x 2 )- N h(x)} 


{(l+x 2 ) N F(x)J- dx. 


Since (l+x 2 )"\(x) is absolutely integrabld and (l+x 2 ) N F(x) is a 
test function, this integral exists in the ordinary sense. 

Definition 7 . If h(x) e I, then the image of h(x) in S is called 
the generalized function defined by h(x) and is denoted by the 
same symbol. 

Thus each ordinary function h(x) e I defines a generalized 
function and the integral (2.19) has meaning in the ordinary sense 
and in the generalized theory, and in this case these meanings are 
the same. 

Using (l.l8), ( 1 . 23 ) and ( 2 . 5 ) > with definition 6 we deduce: 


Theorem 3 . If f(x) e S and Z(f) = g(y), then 

2j(i ? 

(a) Z[f(ax+b) ] = g(|), a =}= 0, 

(b) Z[f'(x)] = 2rtiyg(y), 

and 

(c) Z[g( -x) ] = f(y) . 

The following two theorems eliminate the possibility of con- 
fusion in the notations f(x) and f'(x), where f(x) and f’(x) are 
contained in I. 

Theorem k . If f(x) and f'(x) exist as ordinary functions contained 
in the set I, then the derivative of the generalized function defined 
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by f(x) is the generalized function defined by f'(x). 


Proof : In the generalized theory, with f(x) and f'(x) interpreted 

as generalized functions, we have shown that 


f f’(x)F(x)dx = -f f(x)F’(x)dx, (2.20) 

-oo -OO 

for all test functions F(x). 

Considering f(x) and f'(x) as ordinary functions the integrals 

oo 


oo 


/ f'(x)F(x)dx, 
-oo 

oo 

J f(x)F(x)dx, 
-00 


and 


oo 

/ f(x)F'(x)dx 
-oo 

each exist. Integrating the first integral by parts we have 

r°° b r°° 

/ f’(x)F(x)dx = lim f(x)F(x) | -/ f(x)F’(x)dx. 

'■'-oo fT°° -a -oo 


b->oo 


Hence 


and 


Ji# -oo 


J3j6oo 


OO 

must be finite. But the existence of/ f(x)F(x)dx implies that both 

-oo 

limits are zero. Hence (2.20) holds for ordinary functions f(x) and 
f'(x). 


Theorem 5 » If f(x) is em ordinary function which is absolutely 
integral on (-oo, oo) — so that its Fourier transform g(y) exists 


16 


loy the classical Fourier integral theorem— then the Fourier transform 
of the generalized function f(x) is the generalized function g(y). 

Proof; We have 


p -1 -1 00 

/ |[l+y 2 ] g(y)|dy= / |[l+y 2 ] f f(x)e" 2jtixy dx|dy 

-00 J -00 J -oo 


00 00 

< [/ ti+y 2 ] dy] f | f(x) | dx 
—00 J-oo 


< 00 . 

Hence g(y) satisfies definition 7 ? and the generalized function 
g(y) exists. Now both the ordinary function g(y) and the generalized 
function g(y) satisfy 


00 00 

/ g(y)G(y)dy = / f(x)F(-x)dx 
-co —00 


for any test function F(x), where Z[F] = G(y). Therefore , the 
Fourier transform of the generalized function f(x) is the one 
defined by g(y). 


2.2 The Dirac delta function 




1. The sequence |(^) l "e represents the important Dirac 
delta function S (x) which has the property that for any test 
function F(x), 

CD 


/ 8(x)F(x)dx = F(0). 
J-00 

To show this, we note that 

n™ i 2 

(a) J (S) 2 e- ra 


( 2 . 21 ) 


dx = 1; n = 1, 2, . . 


-00 
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00 1 


dx = (i-) 2 . 


(») f 0 2 M 

-00 

If F(x) •€ S, using (a) we write 
oo 1 2 

| J (|) 2 e _nX F(x)dx-F(0) 


rot 


oo 1 


/ (^) 2 e" nX [F(x)-F(0)]dx 
-00 


oo 1 


r (*) 2 e -™ [F(x)-F(p ll xdx 
J 'n x 


-oo 
oo 1 


</ (2) 2 |x|. 


-nx 


-oo 


max |f’(x) 
- 1 
(«n)2 


F(x)-F(0l 


0 as n — > oo 


where in the last line we used the Mean Value theorem for derivatives, 
i.e., on any interval [0,x] there exists an a, 0 < a < x, such that 

F .( a ) = F(x)-F(p _ ) _ and thus | F '(x)| > |F’(a)| on any interval 
[0,x]; and (h). 


2. The Fourier transform of $(x) is 1. To see this we note that 

1 


2 2 

-* y 


~ 2 ^ n 

,m2 -nx < > e , 

V e 


and hence the sequence 


m 


is a representative of the Fourier 


transform of $(x). But for any test function G(y) 

2 2 _ 

-« j- r °° 

45 oo J e n G(y)dy = J l*G(y)dy, 


-oo 


-oo 


hence 


$ (x) < > 1. 
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-2jtix y 

3. The Fourier transform of $(x -x q ) is e 0 . This follows 
easily from part (a) of theorem 3. 

-2«ix y 

4. From $(x-x q ) 4- > e and part (c) of theorem 3> 

we have 

2«ixy 

e 0 4 » S(y-y Q )- 

For y Q = 0, we have 

1 < > $(y)* 


^ 2rtixy -2jtixy 

5 . Writing cos2«xy Q = ^ [e °+e °] and applying the 

above relation gives 


cos2«xy Q < >|[ $(y-y Q ) + 8(y+y 0 )l* 


. 2itixy -2jtixy 
Writing sin2nxy = rr [e u 0 


'o 2i 

same relationship, we find that 


- e 


] and using the 


sin2nxy < > rr ( 8 (y-yJ - fc(y+y.)]' 

O o o 


Definition 8 . Let g(x) be an ordinary function such that, for any 
test function F(x) which is zero outside of (a,b), g(x)F(x) is inte- 
grable on (a,b), a < b. If f(x) is a generalized function such 
that 


co 00 

/ f(x)F(x)dx = / g(x)F(x)dx 
-00 J -00 


then we define 


( 2 . 22 ) 


f(x) = g(x), 

for a < x < b. 


In the sense of this definition we have 8 (x) = 0 for 0 < x < 00 
and -00 < x < 0. Suppose F(x) = 0 when 0 < x < co and when -00 < x < 0, 
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then by continuity F(O) = 0 and 



dx = 0, 


2.3 The Convolution theorem . 

In section 2.1, we have shown that (1.17) - (1*22) hold for 
generalized functions. It can he shown that (1.23) and (1.24) also 
hold. However, in order to extend (1.26) and (1.28) to include 
generalized functions, some restrictions must he placed on these 
functions. If we try to give the symbol 

00 

(f * g)(x) = f f(z)g(x-z)dz ( 2 . 23 ) 

-00 

a meaning for generalized functions f(x) and g(x), then in the case 
where f(x) and g(x) are ordinary functions, ( 2 . 23 ) must have the same 
meaning as ordinary convolution. We have already found that for the 
convolution theorem restrictions are necessary on the ordinary functions 
f(x) and g(x). 

Suppose the generalized functions f(x) and g(x) are defined 
by the sequences { f n ( x )} and {s n ( x )} respectively. Consider 
the sequence 


a (x) = 
n 



x-z)dz. 


(2.24) 


For each n, a (x) is a test function. To see this, we let f (x) < > 

7 n n 

F (y) and g (x) < > G (y). Applying the classical convolution 

n n n 

theorem gives 


00 

a (x) = f F (y)G (y)e 2ltixy dy. 
n J -oo n n 
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Now P (y)G (y) is a test function and, by theorem 1, so is a (x). 
n n n 

For H(x) e S with Z’^H) = h(y) 

r°° r 00 r°° 

J a n (x)H(x)dx = J H(x)dx J F n (y)G n (y)e 2nixy dy 

-00 -00 -00 

00 00 
" J F n (y)& n (y)dy J H(x)e 25tixy dx 
-00 -oo 


-f F n (y)G n (y)h(y)dy. 
-00 


Now h(y) is a test function, but in general the sequence ^F n (y)G n (y)|. 
is not regular. Hence { a n ( x )} W not define a generalized function. 


In the functional analysis approach (see [4], p. 106 ), it is 
shown that by appropriately restricting f(x) ,(f*g)(x) is a generalized 
function for all g(x) e S. Of primary interest to us here is the case 
when f (x) is a linear combination of Dirac delta functions. To avoid 
some tedious convergence problems which arise in a general approach, 
we shall restrict our discussion to this case. 

x)^ and 

define S(x) * g(x) by the sequence 


Let g(x) be any generalized function defined by K ( 


a (x) = / 8(z)g(x-z)dz. 
n -oo n 

Since g n (x-z) is a test function, we have, by (2.21), 

a n (x) = gjx). 


Also 


oo 

a n (x) = J g n (z) 8(x-z)dz, 


-00 
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Hence 


fc(x) * g(x) = g(x) * 8(x) = g(x), for all g(x) e S. 
Since Z( 8 ) = 1, we have 

Z( I (x) * g(x))= Z(8)Z(g). 

In a like manner we find that 


Since 


and also 


we have 


$(x-x q ) * g(x) = g(x-x o ). 

-2nix y 

Z(g(x-x Q )) = e ° Z(g) 

-2«ix y 

Z( $(x-x q )) = e ° , 

Z(>(x-x q ) * g(x) ) = Z( 8 (x“X q ) )z( g(x) ) . (2.25) 



Proof ? This follows easily from (2.25), and the linearity property 
of the convolution and of Z. For we obviously have 

(*! SU-x^+ag 8 (x-x 2 )] * g(x) = a^ ^(x-x.^ * g(x)]+a 2 [ 8(x-x 2 ) * g(x)]. 
By induction 


N 

A(x) * g(x) a^. [ S(x-x^) * g(x)]. 

j=-M 


22 


Applying Z to both sides of the above equation and using its linearity 
property gives (2.27). 

Let H(y) = Z(a(x)) and G(y) = Z(g(x)). Noting that 

N 

A(x) * g(x) =y a^g(x-x^) 
j=-M 

and applying Z 1 to both sides of (2.27), we have 
N 

y a.. g(x-x^) = Z" 1 (H(y)G(y)). (2.28) 

j=-M 


2.4 Trigonometric series . 

If f (x) is a generalized function for each value of the 
parameter z and if f(x) is a generalized function such that 

oo oo 


.im / f (x)F(x)dx = / f(x)F(x)dx 
^vLoo z -oo 


(2.29) 


for all F(x) € S, the f (x) is said to converge (weakly) to f(x) and 
we write 

lim f (x) = f(x) . 
z->a z v v 

With this definition of convergence in S, we have the following 
theorem. 


Theorem 7 . 


The trigonometric series 


oo 


1 


. X 

in*- 

P 


a e 
n 


n=-oo 


(2.30) 


converges in the sense of (2.29) to a generalized function f(x) if 
and only if a R = 0(|n| N ) for some N as |n| -> oo. If ( 2 - 30 ) converges, 
then its Fourier transform is 
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(2.31) 


00 

g(y) 


n=-oo 

Also f(x) =0 only if a^ = 0 for all n. 

For a proof of this theorem, we refer the reader to [3)> PP* 

58-60. ^ 

The function g(y) is called a "row of deltas" of spacing 
This comes from the equality on any interval (§p>|~) of s(y) 

and an ordinary function which is zero on (^p’l^)* If f ( x ) is 

a periodic and has a Fourier series representation, then the a n 
are the Fourier coefficient* 

p 

p -in*x 

a = 5- / f(x)e Pdx. 

n 2pJ p 

This is equivalent to the statement that convergence in the ordinary 
(strong) sense implies convergence in the sense of (2.29) and that 
the l imi ts are the same. The converse is not true, for by theorem 7 

00 


I 


cos rot — converges in the sense of (2.29) > but obviously not to an 
P 

n=-oo 


ordinary function. 

2.5 The time-frequency interpretation of x and y. 
If we take x to be time t, the function 

2jtif t 

g(t) = e 0 

of frequency f > 0 has transform 


G(y) = $(y-f Q ). 

In the sense of definition 8, G(y) = 0 for y =j= f Q , i.e., G(y) displays 


2b 


the frequency of g(t). If 


-2* if t 

I(t) 

then 

G(y) = $(y+f o )> 

and G(y) = 0 for y ^ -f Q * Hence, if we admit negative frequencies, 

G(y) displays the frequency of g(t) , and we are led to interpret y 
as frequency. 

If h(t) is a function which can be written as 

2nif t 
n 

n=-M 

where M,N may be infinite provided the series converges in the sense 
of (2.2 9), then 

N 

H(y) V(y-f n ). 

n=-M 

H(y) = 0 for y + f , -M < n < N. Thus H(y) displays the frequencies 

n * — 

f and the complex amplitudes a^ of the components of h(t). 

A Fourier series has for its transform a "row of deltas" 

(2.3l), and thinking of the Fourier integral as the limit of a Fourier 
series as the period 2p tends to infinity (see [8] or [9 ] ) > the 
frequency interpretation of y carries over into integral sums. Using 
the symbol f instead of y, we have that if h(t) is a function of time 
with Fourier transform H(f), then H(f) displays the frequencies and 
the complex amplitudes of the components of h(t). 


N 

h(t) a n e 
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CHAPTER III 

\ 

FILTERS 


3.1 Linear systems . 

A linear system, for our purposes , is a linear operator L 
which maps S into S. By linear we mean that for all f(t),g(t) 
e S and all scalars a,b 


L[af(t) + bg(t) ] = aL[f(t)]+ bL[g(t)]. (3.1) 

Let L[ $(t-e)] = h(t;e), and suppose that the family of functions 
h(t;€) is known. Let g(t) be an arbitrary function which we will 
refer to as the input to the linear system L, and let f(t) = L[g(t)]. 
f(t) is called the output of L. Now since 

00 


g(t) = / g(e) $(t-€)de, 
-00 


we have 


f(t) = L [ J g(e) g(t-e)de]. 

-OO 

Assuming that ( 3-l) is sufficient to write 

OO 00 

L[ J g(e) 8(t-e)d€] = J L[g(e) S(t-€)]d€, 


-00 


-00 


then 


Suppose that L satisfies the condition: 

(A) If L[g(t)] = f(t) and t Q is real constant, then L[g(t-t Q )] 

= f(t-t ), i.e. , L is t ime - invariant . 
o 

Then if L[ *(t)] = h(t), L = h(t-t o ) and 


oo 

f(t) = / g(e)h(t-e)de, 
J -co 


(3-2) 


that is, the output of L is given in terms of the input and a unique 
function h(t). The function h(t) is called the impulse response or 
weight function of the linear system L, and its Fourier transform 


H(f ) = f h(t)e" 2 * ift dt 


(3-3) 


H.00 


is called the system or transfer function of L. 

We note here the different uses of the symbol f . When used 
alone or as an argument, f denotes frequency. (See 3 . 3 ). When f 
is written with an argument, f(t), it denotes a function of the time 
variable t. (See 3 . 2 ) 

Note that (3*2) is the convolution g*h. If f(t) < — > F(f) and 

g(t) < > G( f) , then using ( 3 . 2 ) and assuming that the convolution 

theorem holds, we have 


f(t) = g(t)*h(t) G(f)H(f) 

and 

F(f) = G(f)H(f), (3*4) 

OO 

f(t) = f G(f)H(f)e 2 nirt df. 

-00 

That is, the Fourier transform of the output of the linear system L 
is equal to the product of the transforms of the input and the weight 
function h(t). We also note that if G(f) is the transform of an input 
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and F(f) is the transform of a desired output, then from (3.4) the 
transfer function of the linear system L giving the desired output is 


H(f) (3-5) 

H( f ) may in general be complex (see (1.4)] 

H(f) = A(f) e 10 ^ 

where A(f) and ©(f) have already been defined in the classical case 
as the Fourier spectrum and phase angle of, h(t), respectively. 

Definition 9 » A linear system L which satisfies (A) is called a 
filter if A(f) is small in same sense on certain parts of the frequency 
axis. A low-pass filter is a filter for which A(f) is small for |f| > 
f where f^ is called the cut-off frequency . A band-pass filter is a 

filter for which A(f) is small outside the intervals [-f ,-f ] and 

c c 

[f c ,f c ]. A frequency f is said to be passed by a filter if A(f) is 
not small. 

3.2 Ideal low -pass filters . 

We shall restrict our attention here to some particular cases 
where the phase single ©(f) is constant, 

©(f) = a. (3.6) 


Ideal smoothing filter . 


This, by definition, is a low pass filter which passes all 
frequencies f such that |f| < f without change and deletes all 


frequencies greater than f c » No phase shift is involved, and 


hence a = 0. 
Thus 


H(f) = A(f) = 


See figure 3 . 1 . 


1 

0 


M < f c 

| f | >f c . 


(3-7) 
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The corresponding weight function is 



sin 2«f t 
c 

jit * 

or in terms of the angular frequency w = 2nf, 


(3-8) 


b(t) = 


sin w t 
c_ 

rtt 


(3 -9) 


If g(t) is the input to this filter, then the output is 


00 

f(t) = / g(z)h(t-z)dz, 
-00 


which has transform [see (3.4)] 

|G(f)H(f) 

F(f) = 


l f l < f c 


f| > 


where g(t) < > G(f). Hence f(t) contains no frequencies greater 


than f . 
c 
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Ideal smoothing and differentiating filter . 

By theorem 3(h), if g(t) < — > G(f), then 


g'(t) <■ — > (2*if) G(f). 
Repeated application of this theorem yields 


g (n) (t) < > ( 2 «if) n G(f) (3-10) 

til 

From (3«5) we see that to find the n — derivative of an input g(t) 
the transfer function must he (2«if) n . Then, in order to smooth 


using the ideal filter and find the n — derivative, the transfer 
function is given hy 


H(f) = 


f(2*if) n 


0 


l f l < f c 

Ifl > f c 


(3.U) 


and the weight function is 


f 

r c 

h n (t) = J ( 2 jtif) n e 2 ]tlft df. ( 3 - 12 ) 

-f 

c 


But differentiating ( 3 . 8 ) n times, we have 


h (n) (t) .J ( 2 «if) n e 2 ” ift df 


-f 


and so 


h n (t> = h (n, ( t ) 


(3.13) 


( 3 - 1 ^) 


Thus to find the weight function of the ideal smoothing and different- 
iating filter we simply differentiate the weight function of the smoothing 
filter the appropriate number of times. 
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Then 



(t-z)dz. 


(3-15) 


3.3 The sampling theorem * 

Ideal filters of the type discussed above are not physically 
realizable because of the jump discontinuities at + f c * Futhermore, 
in digital filtering the input consists of a finite number of equally 
spaced values g^, M < m < N, which we may assume are samples of some 

function g(t) for t = mAt = p-. We may also assume that g(t) defines 

s 

a generalized function, for, recalling theorems 2 and 5, this does not 
place a serious restriction on g(t). It is obvious that g(t) is not 
uniquely determined by the values g^, and hence the set of values 

are associated with a subset 

If, we know that the samples g^ arise from a function g(t) 
whose transform G(f) is zero for |f| > f Q , then the subset of 

S is reduced to a subset G^ c. G^. In this case g(t) is said to 

be band-limited . 

Theorem 8 . Shannon* s ser.pling theorem ( see [ 5 ] ) . 

If g(t) Is band-limited, i.e., if g(t) < > G(f) where 

G(f) = 0 l f l£ f a (3«l6) 

then g(t) can be uniquely determined from its values 


SiN 0f S ‘ 


g n = S ^2f J 


(3.17) 


at a sequence of equidistant points of distance pp- apart. 

a 
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Futhermore 


oo 


e(t) 


sin n (2f t-n) 


a 


g n it (2f t-n) 


n=-oo 


Proof : We first compute the g R . We have, using (3.16), 


a 

g(t) = f G(f) 

a 


2rtift 
e or. 


hence 


f f 

a nitij 


= = / G(f) 

- f a 


a 


df. 


Expanding G(f) in a Fourier series on we have 


00 -ruti- 


G(f) G n e 

n=-oo 


a 


> - f a <f< V 


where 


°n ' 4 / G ' f) 

-f 

a 


f .f 

a niti— 

e a df 


comparing (3*19) and (3»2l), we have 


6 

G n = 2f 


a 


The function 


G(f) = ) 


00 -mti— 

V «n . f a 

e , -OO < f < OO, 


2f. 


a 


n=-oo 

is the periodic extension of G(f) and 


(3.18) 


(3-19) 


( 3 - 20 ) 


(3-21) 


Hence we may write 

G(f) = H(f)G(f) 

where 

*1 


H(f) = 


0 


Now [see ( 3 . 7 ) and (3*9)] 


< f c, 


f > f. 


a 


So we have 


sin 2rtf t 

-VF^" < «< f >- 


cd g 


G(f) = H(f) y 


a 


n=-oo 


and 


00 


n=-oo 


g 

^(f)e 

a 


00 00 


.f 

-niti— 


g(t) = f [ y ^_«(f)e a ] 

-00 ^ a 


2« 1 ft df 


n=-co 


od _ 

■ ikl 


.f 

00 -rnti— 

H(f)e a e 2 nift df. 


00 

n=-oo 

Applying the t domain shifting theorem gives 

00 g _ Sin 2 *f a (t- £-) 


g(t> = y 2? 

L- j / 


n=-oo 


a «(t- |j-) 


a 


00 


n--oo 


sin jt(2f a t-n) 
g n «2f Q t-n 


( 3 - 24 ) 
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If f is any number such that f > 2f , then the theorem remains 
s s — a 

true if in the proof the periodic function G(f) is assumed to be of 

f„ f, 

period f , and Il( f ) = 1 for |f j < s— ; H(f) = 0 for jf j > 5^. 


Therefore 


g(t) = ) g. 


sin tt( f t-n) 


n u-n; 

. s 


where 


s n = g( f“)* (3 ’ £ 

s 

If the g are known, as assumed above, for M< n< N, then the 
n — — 


function 


h . . 

wM sin n(f t-n) 

= )_, e n rt(f t-n) 
n=M 

differs from each function g(t) of by 


sin «(f t-n) 


€ g^ Z_, «(f t-nj 

n=-co s 


sin «(f t-n) 


r-> sin 
/ , g n n(f t-n 


n=N+l 


where g^ = g(— ) . Hence, at least in the cases where the series in 
s 

(3.25) converges uniformly to g(t), the maximum difference 


Max |e ^(t)|= max IsUJ-fi^t) | (3-29) 

t t 

can be made as small as we please by taking a sufficient number of 

terms in g^Ct). Hence we can associate with the samples W a 

unique function g(t) in the sense that (3.29) can be made arbitrarily 
small by talcing a sufficient number of samples. 
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3 * J + Definition of a digital filter . 

Suppose that the sampled function g(t) is hand-limited. 

Then G(f) =0 for |f| > f Q . If H(f) is a desired transfer function, 

then H(f)G(f) = 0 for |f| > f^. Thus if H(f) is a periodic extension 

of H( f ) with period f > 2f , we have the transform F(f) of the output 

f(t) given by 

F(f) = H(f)G(f) = H(f)G(f), (3.30) 

for all f. 

If H(f) is such that H(f) can be written as a trigonometric 
series, 

00 0 „ .f 

V -1 2n « 1 

n(f) =2 j v s (3 * 3l) 

n=-oo 

with a n = 0(jn| N ) for some N as |n| ->cd # then, by theorem 7, 

H(f) e S and is the transform of 

00 

h(t) = V a fl 8(t+|- ) (3.32) 

L 1 g 

n=-oo 

Now g(t) is time sampled. In order to obtain a time sampled version 
of the output f(t) we might try to define a convolution h(t)*g(t ) 
and extend theorem 6 to functions A(t) = h(t). Assuming that we 
could do this, we would have 

00 

f(t) =\ a n g(t+|-) 

' s 

n=-oo 

which wo\ild yield the sampled version of f(t) for t = as 
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which is impossible to use digitally since it requires infinitely 
many samples. 

Alternately, let 

N 2mtij- 

W f) ■ I v s (3 - 33) 

n=M 

be a trigonometric polynomial which approximates H(f) in some sense. 
Then ( 3* 33) is the transform of 

N 

= X a n (3-3^) 

n=M S 

and we have shown that the convolution h^Ct) * g(t) is defined for all 

g(t) e S. Also, theorem 6 holds. Thus 

F(f) = G(f)H(f) = GCfjH^Cf) = F(f) (3.35) 

and f(t) < — > F( f ) is given by 


f(t> = 


N 


g(z) S( t-z + ~) dz 


n=M 


N 


a n g(t+ fT )- 


n=M 

For t = 2-, f = f( S-), g = g( |-) , we have 
f s m f c m f s 


N 


f = 
m 


i n®mfn 


n=M 


(3.36) 


diis is the fundamental formula of digital filtering. 


I 
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Note that any pair (3*33) and. (3-34) determine a linear operator 
L on S which satisfies condition (A) and, on the subspace G g of all 
band-limited functions g(t) with 2f < f , acts as a low -pass filter. 

Now any finite set of constants a determines a function ( 3 . 34) , 

n 

which determines ( 3» 33) and hence a linear operator L. 

Definition 10 . Let a , , M < n < N, be any set of constants. Then 

the linear system L determined by the a is called a digital or 
numerical filter . 

Application of L must be limited to the subspace G g . Otherwise 

"frequency folding" occurs, i.e., frequencies in the intervals 
(2n-l)f (2n+l)f 

( g , 2 n = ± 1, + 2, . . . are folded back into the 

f f 

( !■, - 77 ) • For example, suppose the input contains a frequency 

f 

component Acos 2 jt( f^+kf^Jt where f q < -| and k is a positive 
integer. Then if we sample at t = 

s 

Acos [ 2 «(f +kf )^-] = Acos [2jtf + 2 nkn] 
o s x 01 

s s 

= Acos (2jt f Q ~) . 

s 

The sample values would be the same as those obtained from a component 
Acos 2«f t for t = . Hence the filter treats the frequency 

OI 

s 

f 

s 

f + kf > — in the same manner as f . 
o s 2 o 

3.5 Even and odd transfer functions 

In cases of interest here, the transfer function H(f) is either 
even or odd. Hence the trigonometric polynomial H^f) which 

f f 

approximates H(f) can be written in terms of cos 2mt— and sin 2rot— 

s s 
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respectively. If we take M = -N some advantages are gained. Let 


N 2nnl-, 


•w f) = v f) 


„ v 


/ a e 

L n 

n=-N 


For even functions, 


N 


V f) ■ % + 2 I 


a + 2 ) a cos 2rut— . 

° 4 1 n f s 

n=l 


For odd functions. 


N 

V f ) = -I 


2i ) a sin 2rut— 
n i 

n=l £ 


uoi; 


( 3 - 38 ) 


( 3 - 39 ) 


Two questions now arise: 

(1) given H(f), how are the coefficients a n to he chosen, and 

(2) what is the error introduced by the approximation 
F(f) = F( f ) ? 


3.6 Methods of filter approximation 

If H(f) is an ordinary function, there are currently two methods 
of approximating H(f) and obtaining the coefficients a n » One of these 
methods — the Min-Max technique — is given by Martin [ 7 ]. Essentially, 
it assumes continuity of H(f) in which case, if Q n (f) as a se ^ N 

f f 

continuous and linearly independent functions on [ then 

there exists a polynomial 


Pjj(f) ~ + . . . + a^Q^(f) 

f f 

which deviates the least from H(f) on ( 1-, ~) , i.e.. 
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N 


max 

f 

f e (- 


f 

s S\ 
2 ,— 2 ' 


|H(f)- 


■p N ( f) l 


max 


f € (- 


f f 

S S\ 

2 » 2' 


|H(f) 


x q (f) 

n n 


n=l 


for any numbers x i» * 2 * • • ^ 1116 Q^(f) are obtained after 

putting a constraint (or constraints) on a trigonometric polynomial 
(3-37). P N (f) is then fitted at a finite number of points to H(f) 


in the above sense. A good approximation of the a is obtained by an 

n 

iterative process, but the technique is long and complex, and not 

very versatile. That is, any change in H(f) necessitates a complete 

repetition of the process for finding the a . 

n 

The alternate method 'assumes that H(f) can be approx imat ed by 
a Fourier series. 


where 



(3-^0) 


(3-41) 


and is taken to be the truncated series for H(f), 

N 2mri|- 

V f) - X h n e S ' (3-te) 

n=-N 

This gives a function which is the best fit to H(f) in the least mean 
square sense. 

Noting that, since H(f) =0 for |f | > 
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(3-43) 


h(t) =J H(f)< 


J. f 

2 2ruti— 


df, 


- f 


_s 

2 


and comparing (3*43) and (3»4l), we see that 


“„ = F *-£■>• 

S S 

This is the basic formula for computing the h = 


n 


(3-44) 


a to use in ( 3« 36 ) * 
n 


Therefore (3.36) can be written as 

N 


f m X ^n^m+n . 


n=-N 


(3-45) 


3.7 Error analysis . 

With an approximation Hjj(f) of H(f), (3*35) becomes 


F(f) = G(f)H(f) = G(f)H N (f) = F(f), 


and so 


F(f)-F(f) = G(f)[H(f)-H w (f)]. 

This gives the error between the spectrum of the desired output and 
the spectrum of the actual output. 


2* if t 


For a complex frequency component g Q (t) = Ae 


in the input 


we have 


2nif t 

g (t) = Ae 0 < =► A 8 (f“f ) = G(f ) 


and 


F(f Q ) = A S(f-f Q )H(f o ), 
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also 


p(f o ) = a 

Denoting the difference in the outputs by e(f Q ,t) we have 

oo 

I e(f Q ,t)| = \Ja S(f-f o )[H(f o )-H N (f o )]e 2rtift df 
-00 

2«if t 

= |Ae ° [H(f o ) - H^fjll 
2nif t 

= |Ae ° | . | e(f o ,N)|, 
where e(f Q ,N) = H(f Q ) - H^). 

In the time sampled version: 

2rtif Q j- 

I e(f , §-)| = |Ae s | • | e(f ,N) | . (3-46) 

s 

Thus the magnitude of the error in a component of the actual sampled 
output is given in terms of the magnitude of the corresponding component 
of the input function, and of the magnitude of the error in the approx- 
imation of H(f). 

Approximations of e, 

e = max | e(f,N)| = max |H(f) - HJf)| (3.47) 

f f " 

derived mathematically are usually found to be so large as to render 
them useless in applications. In applications of the filters dis- 
cussed later, acceptable values of e are in the range .005 < e < .01, 
or referred to unity, ^%and 1%. When speaking of percent error 
we will always mean e referred to unity. For a given H(f), an N is 
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found empirically such that H^f) approximates H(f) within the desired 
limits. 

However* satisfying the requirement that .005 < £ < .01 does not 
imply the output error is within these bounds (see Chapter VIl). 

3.8 The Gibbs ' phenomenon . 

When approximating an ideal or designed transfer function 
H(f) having one or more jump discontinuities with a truncated 
Fourier series, there exist oscillations in the approximating 
transfer function H^f) near the discontinuities of H(f) due to the 

Gibbs' phenomenon (see [10]). No matter how large N is taken, e 
cannot be brought within the acceptable range .005 < e < »01. To 
avoid this difficulty H(f) is first approximated by a function which 
is continuous. In most cases, this imposes a restriction on the 
input g(t). The particular cases of interest here shall be dealt 
with in the next chapter. 


CHAPTER IV 
FILTER DESIGN 

4.1 Assumptions about the input . 

In order to apply a digital filter to a set of data k) » 

we have made two assumptions about the data: 

I. It arises from a function g(t) which defines a generalized 
function, and 

II. g(t) is band-limited. 

In many cases of interest, the Fourier spectrum G(f) of a signal 
g(t) consists of a desired signal spectrum in an interval [-f ,f ], 

an unwanted signal spectrum (noise spectrum) in intervals t-f a ,f c ) 

and (f ,f ] , and G(f) =0 for |f| > f Q . When applying a low-pass 

filter, elimination of the unwanted spectrum is desired. Hence the 
ideal filter transfer function, H^(f), is such that H^(f) =0, 

|f| >f • Usually H ( +f ) 4 0 and H (f) has jump discontinuities 
c x *" o • X 

at f = + f . If the truncated Fourier series of H T (f) is used to 
— c i 

approximate H^(f), then, due to the Gibbs' phenomenon, large 
oscillations persist in a neighborhood of + Futhermore, the am- 

plitude of these oscillations remains constant with increasing N. 

The truncated Fourier series is continuous everywhere because it is 
a finite stun of everywhere continuous functions. Since H^(f c ) 4 0, 

we expect that the truncated series, H^( f ) , is such that H^f^) 4 °* 

Then, by continuity, H^( f ) is non-zero on some interval (f c> f c +Af) 

where Af > 0 and depends on N. Any unwanted frequencies which 
appear in this interval are passed — though somewhat attenuated — by 
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the approximating filter. Hence, in addition to the large oscillations 

which appear near +f , unwanted frequencies arbitrarily close to +f 

c c 

cannot be eliminated by increasing N. This undesirable property must 

be tolerated because it is a property of any truncated Fourier series 

such that H^(f c ) ^ 0. Hot/ ever, the large oscillations are caused 

by non-uniform convergence of the Fourier series of H^(f). This can 

be remedied by redefining H^(f) so that it is a continuous function. 

We choose to do this on the intervals [-f c ~Af, -f ) and ( f c> f c +Af) 

for some Af > 0. Any unwanted frequencies in these intervals will 
be passed to some extent by the filter, but, as pointed out above, 
this cannot be avoided anyway. However, in many applications unwanted 
frequencies do not appear near ±f c * Therefore, we make the following 

third assumption about the data: 

III. The desired signal spectrum and the unwanted spectrum of g(t) 
are disjoint. 

Then there exists a Af > 0 such that the signal spectrum G(f) =0 
on ( -f^-Af , -f c ) and (f , f^+Af). Letting f^ - f c +Af , we may modify 

IL(f) on [-f m , -f ) and (f , f ] to obtain a function H(f) contin- 
1 1 c cl 

uous for all f and thereby eliminate the Gibbs' phenomenon. H(f), 
as defined on the intervals [-f T , -f ) and (f , f^], is called the 

roll-off of the filter, and the frequency f is called the termination 
frequency . 

4.2 Filter design by convolution . 

The usual approach to the design of a filter is to select the 

ideal transfer function H (f) on [-f , f ] and then to specify the 

I c c 

roll-off. This gives the filter transfer function H(f) from which 
the weight function h(t) is found. The weights of the filter to 
be used in (3*^5) are then computed from (3*44). In addition to 
not being very versatile, this approach usually involves some 
rather long and tedious integration in determining h(t). 
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We propose a different approach to the design which simplifies 
the integration and gives considerable freedom in varying the roll- 
off shape of the filter. We shall use the convolution theorem of 
Chapter I: 


oo 

g(t)k(t) < > / G(f-z)K(z)dz 

-GO 


(4.1) 


where g(t) < > G(f) and k(t) < > K(f). 

We note here that filters for simultaneously performing smoothing 
and differentiation can be found from the weight-transfer functions, 
h(t) and H(f), of the smoothing filter in a manner analogous to that 
in the ideal case [see section 3.2]. That is, to smooth and find 
"fcht 

the n — derivative, the transfer function is 


With 


we have 


where 


Af) = (2«if) n H(f). 

y n (t) <= > Y®(f), 

y”(t) = h (n) (t) 


h(t) < > H(f). 


In (4.3a), let t = =£. Then t n = and dt n = (- ^"dx". 


(4.2) 


(4.3a) 


Hence 


n, -x^ 

y (7-) 


dV?) 

s 

s 


d"h(^) 


= (-l) n f n — 

' s „ n 


dx 
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Using ( 3*44) to compute the weights of the filter, we have 


v n = L. y * ( ±) 

A f ' -f ' 
s s 


= ^-y n (— ) ] 

r v f x=k 
s s 


= (-l)“f 

We now see that we may write 

n / . >n n d \ 
yv = (- 1 * ) f — zr 


{a n ^(^)\ 

n l s . s Ji 


dx‘ 


n 


^x=k* 


n 


s dk n 


(4.3b) 


X -k 

where h^ = ^r-h(— ) and, for purposes of differentiating, k is treated 
s s 

as a variable in the right side of (4.3b). 

Returning to the problem of designing the filter, we conclude 
from the above that we may restrict ourselves to the design of smooth- 
ing filters. Hence suppose that 

oo 

H(f) =J^G(f-z)K(z)dz. (4.4) 

Ideally, for smoothing we want H(f) to be continuous, and 


H(f) 


1, 0 < f < f , 

7 — “ C* 7 


^monotonic decreasing, f^< f < 


0, f > f c , 
H(-f), f < 0. 


We attempt to find function G(f) and K(f) such that H(f) 
given by (4.4) has these properties. 


(4.5) 
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Then the weight function h(t) is given by 


h(t) = g(t)k(t) . 


Let 


G(f) = 


1 , 


f < 


f 

c T 


- 2 


f > 


f c +f T 


Then comparing with ( 3 . 7 ) and ( 3 . 8 ), we see that 


g(t) = 


sin n(f T +f c )t 


Tit 


Noting that G(f-z) = 0 for |f-z| > — (4.4) becomes 


f T +f c 


f+ 


f m +f 
T c 


H(f) = J K( z)dz. 


f- 


f T +f c 


To find H(f Q ), K(z) is integrated over an interval of length (f, 
with f as its mid-point. Any function K(z) which is zero for 
f m -f 

I I *P fy A+* 

z| > — 5 — = — and is an even function of z with area 1 on 

yield - s a satisfactory H(f). 

Filter 1 . The Ormsby smoothing filter (p=l) . 

In (4.8) let 


(4.6) 


(4.7) 


(4.8) 


■ +f c> 


47 



K(f) = K^f) = 'S 


l^o, |f | > 

See figure 4.1. K^(f) 



Then 



2 _ 

Af 



(2itft)df 


A 

-i& [ sto2 " ft ] 

0 


sin rtAft 



and 


h x (t) - k^tJgCt) 


sin *Aft sin jt(f m +f )t 

T c 

2 2 * 

* Aft 


(4.10) 


Changing to the angular frequency w = 2«f, A w = 2nAf, w T = 2*f T , 


w = 2jtf , we have 
c c 


h x (t) = 


/wt (V^O 

2 sin sin — T . c 


«Awt 


and applying a well-known trigonometric identity 

cos w t - cos w m t 
h(t) = S T. . 

i&rt 

This is the weight function given by Ormsby [ll] for p=l. The 
corresponding transfer function as a function of f is 


(4.11) 


1, 
0, 

f+f 


f -f 
T 

~&T’ 


If I < f c , 

|f| > fm. 


- f T ^ f < - f c> 


f < f < f_. 
c - T 


H^f) has a straight line roll-off (see figure 4.2). 

Replacing f by 2nf = w gives the transfer function as given by 
Ormsby in terms of w. We remind the reader that, when using the angu la r 
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Figure 4.2. 


frequency w = 2«f , a factor of ~ appears in the statement of the 
Fourier integral theorem, i.e., if 


then 



In the statement of the theorem in Chapter I, we have 



G(f) = G(w) = G(2«f). 

For convenience, we shall use the same symbol G for both G and G. 
Then the argument of G determines which form of the theorem to use, 
i.e., G(f) indicates the form of Chapter I is to be used and G(w) 
indicates that the above form is to be used. 

dH (f) 

Finally, note that — — — is discontinuous at + f c and + f^. 
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Filter 2. The Martin-Graham smoothing filter. 


In (4.8) let 


K(f) = Kg(f) = 


Tt Jtf 

2Af Af ’ 


kl <%, 


kl >^§ • 


See figure 4.3» 


K 2 (f) 



Then 


Figure 4 . 3 . 




= - f 

Af J 


cos COS 2rtft df 


sin -2*t)f sinter; +2*t)f 


2(^ - act) 


+2*t) 
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Then 



sin(| - «Aft) 


( 


Af 


2t) 


sin(^ + nAft) 

— T 1 

<S? + 2t) 


1 r cos rtAft 

W-Bt) 


cos nAft ^ 

<s + 2t > 


cos nAft 
(l-4Af 2 t 2 ) 


h 2 (t) = k 2 (t)g(t) 


cos itAft sin «(f T +f c )t 
«t(l-4Af 2 t 2 ) 


(4.13) 


2 2 2 2 
where Af = (Af) . We will also use the notation Aw = (da/) . 

Letting w = 2nf gives 


h 2 (t) = 


Art , <V*c )t 

cos -g- sin 2 

___ 

«d- 

« 


and using a well-known trigonometric identity gives, after simplifying. 


h 2 (t) 


*(sin w c t + sin w^t) 
2t(it 2 - Aw 2 t 2 ) 


(4.14) 


This is the form of the weight function given by Graham [12]. 

The form given by Martin [6], [7] is obtained from (4.13) by 

f Af* 

going to the frequency ratio r = — , 2h = r d = — , etc., and computing 

s s 
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VrVl' 

s s 


cos [*(Arf g )(- —)) sin [«f g (r c +r T )(- 
I s s 


«<- rHi-'fr 2 f 2 a,) 
s f 

s 


cos [njtr ] sin [rut(2r +rj] 
J. c g 

nn(l-4r ,^n^) 

a 


(4.15) 


This is a convenient expression for confuting the weights h of the 
filter. The value of h Q is computed by using L' Hospital' s rule, and 


h = 2r +r 
o c d 


v\ 

* 1 1 « 

2* 


(4.16) 


The same procedure must be used for finding h if m = giving 

m 2r 


d 

o rt(2r>+r) 

Ein( 

d 

The transfer function of this filter, in terms of w, is 

1 M <w p , 

o; 


. n jt(w-w ) 

H 5 (w) = -j^l+cos —], 


w| > w T , 


Jt(w+w c ) 

|[l+cos ], 


See figure 4.4. 

Alternate expressions for the roll-off are 


w c < w < w T , 


-w m < w < -w 
T c 


n(w+w ) 


2 tl+oos ar 


-] = cos (- 


fl(w-w ) 


2£w 


). 


(4.17) 


(4.18) 
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rt(w+w ) „ jt(w+w ) 

|[ 1+ COS —£-l = CO, 2 (-gij. 

d 2 H 2 (w) 

Note that H (w) has one continuous derivative, and ?> — is dis “ 

2 dw 


continuous at + w^, and + 
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Since K 1>lf (z) is the truncated series of a function with jump dis- 
continuities at + — [see(4.9)], the Gibbs 1 phenomenon is present. 

A p 

Hence overshoot is present near + — , the amplitude of which can not 

he reduced by increasing N. We can expect some relatively large 

oscillations to be present, at least for small values of N, in 

%,„«• K 2 (z) is continuous, and the amplitude of the oscillations 

° f K 2’l/ Z ) decreases monotonically with increasing N. Hence we expect 

the Martin-Graham filter to perform better than the Ormsby (p=l) 
filter. The results of comparative programs where the truncated 
series (4.19) and (4.20) were computed at equidistant points indicate 
that this conclusion is true. For e = .01, over 50% more weights 
were required by the Ormsby filter. 

4.4 Some new smoothing filters 

We shall give, without performing the details of integration, 
several new designs which are of some interest. The transfer function- 
weight function pair will be given in terms of the angular frequency 
w = 2itf. 


Filter 3 . Let 


Then 


K^(w) 


= 4 


4rt 2,jrw\ 

s cos W 


0 


M < % 


w > 0 


(4.21) 




2 2 sin f t 


k (t) = p „ — ■ 

3 4it^-Aw 2 t 2 Awt 


h 3 (t) = k 3 (t)g(t) 
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(4.22) 


.2 2 sin ^ sin ^ t 

= ( »%o )( S-T 2 ) 


4* 2 -Aw 2 t 2 


Awitt 


4jt c 


), 2 A 2* 2 1 

4it -Aw t 


h,(t) 


where h^t) is the Ormshy weight function (4.1l). The roll-off of 
H^(w) is given by 


± 2(w-w c ) w T -w 


2it 


sin 


Aw + Aw * 


W < w < w , 


and H^(w) has two continuous derivatives, (see figure 4.5) 


Filter 4. Let 


Kj^(w) =* 


COS^ (2t) > 
2Aw Aw 


0 


w <%, 


. Aw 
w > - 


(4.23) 


Then 


9* 


2 Awt 
it cos -g— 


* 1 ^ “ 2 A 2.2 * / 2 . 2 2 \ * 
9« -Aw t (it -Aw t ) 


and 


h^(t) = k^(t)g(t) 


9, 2 ,” 2 

i o o o n r> i J 


2 2 2 2 2 2 
9« -AW t^ it -AW t 


itt 


9« 


. 2 A 2. 2 2 

9 it -Aw t 


Ut) , 


(4.24) 
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is the Martin-Graham \ 


(w) has three continuoi 


o (w-w ) 

: )/ w ) = JE COS 


W T . 


Aw 

This is shown i_. 



Filter 5. Let 


K ? (w) 

This gives a weight fu 

h (t) = — r[sin ( -— — )t]» [2 sin ~ t - tAw cos^ t] . (4.26) 

5 rtZV^t 

The roll-off of H^(w) , w q < w < w T , is a third degree polynomial and is 
essentially the same as that of H^Cw). 

Using the quantity e defined by (3*46) as a measure of the 
performance of a filter to compare the above filters, one is led 
to the following conclusions: 

1) The Martin-Graham filter gives e = .01 with smaller N than 
any of the others. In fact, out of numerous designs none 
has been found which gives € = .01 for smaller N than this 
filter. The performance of filter 5 is essentially the same, 
the e values differing slightly in the third decimal place. 

2) Filters 3 and 4 give values of e < .005 for smaller N than 
the Martin-Graham filter and filter 5* 

3) In no case did the Ormsby filter perform as well as the 
other filters. 

In comparison with the Martin-Graham filter, the only advantage 
filter 5 has is that no special evaluation for h n » n =}= 0, is required; 

h is the same for all the above filters. In addition to the improved 
o 

performance for e < .005, useable error bounds can be found for filters 
3 and 4 without resorting to empirical methods. 


3* (1 _ 

Aw U 

Aw 


0 , 


I | ^ Aw 

M < -£ , 


W| > ~2 


(4.25) 


3 3 

inction, where Aw = (Aw) , 
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4.5 Some smoothing error bounds 

Except for filter 5> each of the above weight functions are of 
the form 


h(t) = 


pTty 


where k(t) is an expression containing stuns and products of trigono- 
metric functions of t and P(t) is a polynomial in t. The Fourier 
coefficients of H(f) computed from h(t) retain this character, 

k(^) 

h „ = r h( r ) = r — =|- • 

n f s f s f s P(^) 

s 

Now the error as a function of f and N is 
e(f,N) = H(f) - H„(f) 


= 2 ) h cos 2nn — 

L n f c 

n=N+l 


oo k(- £-) 

T Y — r r cos 2nit f“ * 

jL *- F> 

s 


n p 

Letting A = ma^ |k(- — ) cos 2nn — | , we have 


e = max |e(f,N)|< — ^ I^tI- (4.2 7) 

s n=N+l n f ' 
s 


If |p(t) I > 0 for t > the sum in (4.24) can be approximated by 
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Martin-Graham bound 


dx 




The above method gives (see figure 4.6) 


2 2 

^ 1 . 4N Af 
e < - log 


* ~ 4N 2 Af 2 -f 2 


For e = .01, the predicted value of N is 


N > 


2.85 fa 


Af 


and for e = .005 


4f 

N > — - 
- Af 


These values of N are much too large. It has been determined 

1.25 f„ 

empirically that N > — ^ — - gives .005 < e < .01. 

Filter 3 « 

For this filter. 


c ^ 2 


1 

7 

c 

n 


log [ 


NAf+f 

C 

Nftf-f 


-) - 


2f 

s 

NAf 


> . 


For e = .01, the predicted value of N is 

2f 

“ - Af’ 

3f 

and for N > — e < .003 (see figure 4.6). 


(4.28) 


(4.29) 

(4.30) 


(4.31) 


(4-32) 
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NAf /f 


6l 


Filter 4. 


For this filter 

e < gjj- {9 log [4rt 2 M 2 Af 2 -n 2 f s 2 ] - 16 log [2«NAf] 

- log [Un^Af 2 - 9« 2 f 2 l) • (4-33) 

s J 

For e = .01, the predicted value of N is the same as in (4.32). 

3f 

For N > (4.33) gives e < .0014 (see figure 4.6). 

4.6 Smoothing filter constraints 

In general, a signal g(t) may have a polynomial content, and in 
such cases g(t) is not band-limited. Denoting the polynomial content 
of g(t) by P(t), if 

g(t) = i(t) + p(t) (4-34) 

where g(t) is a band-limited function, then the weights can be con- 
strained so that the sampled values P(mAt),At = are passed 

1 s 

without error. 

We recall that the output of a digital filter is given by 

N 

*m ^n^ra+n. 

n=-N 

Applying this to the sampled version of (4.34) gives 

N_ 

= I h n'Vn + P[(mfn)At11 

n=-N 

N 

= X "nVn + X h n P[(wn)Atl - (ll - 35) 

n=-N n=-N 
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Since g(t) is band-limited, the first term on the right side of 
(4*35) poses no problems. We want the second term to be P(mAt). 
Assuming that P(t) is of degree p. 


P 

P(t) = ^ a^. 

We want 


(4.36) 


P(nAt) = ^ a..(nAt)^ 

J=0 

N p 

V~i 

= ) h ) a.(m+n)^At^. 

L n /. V 

n=-N j=G 


Interchanging the summation gives 


N 

^ h n (mfn) J At J . (4.37) 

j=0 j=0 n=-N 

V, T e see from (4.37) that it suffices to consider the k— term 


a.(mAt)'- ) = 


P^ 


a. 

0 


N 

k ,k \ ’ , , *k,k 

m At = ) h n (m+n) At 

n=-N 

k 

or, dividing by At , 


N 

m k = ^ h n (m+n) k . 

n=-N 
k. 

Expanding (m+n) and siunming each term gives 
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N 


N 


k k 
m = m 


^ + krn^" 1 nh fi + . . . 


n=-N 


n=-N 


N 


N 


n=-N 

From (4.38) we see that it suffices to have 
A: 


+( k )m k " r ^n r h n + . . . (4.38) 


n=-N 


N 


h = 1 
n 


(4.39) 


n=-N 


B: 


N 


n J h n =0, j = 1,2, 


>P* 


(4.4o) 


n=-N 


The transfer function of a digital smoothing filter which approxi- 
mates smoothing filters of the types discussed in section 4.2 is an 
even function of f and can be written in the form 


iL. 

H (f) = h + 2 \ h cos 

n o L n f 

n=l 


(4.41) 


The weights are related by h =h . Hence for odd integers j, 

n - n 


n \ = - ( _n)Jh _ n 
or 

n^h + (-n)^h = 0 

n -n 
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and 


N 

n^h n = 0. (4.42) 

n=-N 

Thus (4.4o) is satisfied for all odd integers 3 without imposing any 
conditions on the h n * If (4.39) is satisfied, the filter passes a 
first degree polynomial exactly. If, in addition, (4.4o) is satis- 
fied for o=2, the filter passes a third degree polynomial exactly, 
etc. Practical considerations usually limit 0 to 2, i.e., p=3. 

The simplest way to satisfy (4.39) is to use new weights 



n=-N 


(4.43) 


If N is chosen so that .005 < € < .01, the new weights usually do not 
change e significantly. 

For 0 > 2 the usual approach is to derive the constrained 
weights h n so that the mean square error between the unconstrained 
transfer function H^f) and the constrained transfer function 

N 

H^f) = h Q + 2 ^ h^cos 2mt (4.44) 

n=0 S 


is minimized. 

Note that (4.39) is equivalent to the condition 

V 0) “ 

and (4.40) is equivalent to the conditions 


=0, 1 < j < p. 

f=0 


d%(f) 

df 5 


(4.45) 


(4.46) 
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Taking the case p=3 and using a Lagrangian multiplier, we wish to 

find weights h in terms of the h such that 
^ n n 

f 

s 

N 

\ 2 — 

df +\) n h n 
n=l 

and such that H^(f) satisfies 

aii(f) 

— df + X m 2 . (4.47) 

dh 

m 

The condition (4.45) is incorporated in the following way: 

N 

V°> ■ 5 0 + 2 1 E n " X > 

n=l 

so 

N 

h =1-2 ) h . (4.48) 

o /_, n 

n=l 

Hence 

H^fl-H^f) =1+2 j h n (cos 2n« — - 1) -H N (f). 

n=2 E 


R = / [^(f) - H N (f)]' 


3r 

is minimized, i.e., — — = 0, 0 < m < N, 

dh 

m 

(4.45) and (4.46) for p=3. 

f 

s 

ff- = 2 / tS N< f) - V f>) 


m 0 


Therefore 


and 


m 0 n=l 


cos 2 rut- l) 

s 


N 

h - 2 ) h cos 2rut^-] 

° L n f s 

n=l s 


f 2 

[cos 2nut- l] df + Xm . 

s 


f f s de 

Let 0 = 2n — , then df = " and 

s 


dR 2f s 


dh 

m 


« N N 

J [ 1 - 2^Th n - h^+ 2^ (h^-h^Jcos n0][cos m0-l]d0+ Xm^ 
0 n=l n=l 


2 f * N 

= — - /{[l - 2 ) h - h ] [ cos rae-l] + 2 ) (h -h ) cos n0 
it J [__ i n o , n n 

0 n=l n=l 


cos m0 


N 


-h^) cos n©J- d0 + X nf 

n=l 


rtx , N 

2f . 'r n _ _ . p 

= — §■ \ -n[l - 2 ) h -h ] + * (h -h )) + X m . 
n l n o m m j 

n=l 


Setting this equal to zero gives 


N 


2f T[h - 1 + 2 ) h +h -h ]f+ X m 2 = 0. 
s'- o n m m J 


n=l 


N 


From (4.48) we see that we can replace 2^h n ^ -h 

n=i 
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Thus 


2f jfh -h + h -h ] f + Xm 2 = 0. 
s'. 1 o o m m J 


Let S= h -h , then 
o o’ 


h -h = B 
mm 


X m 

2f * 
s 


( 4 . 49 ) 


Summing both sides of (4.49) from 1 to N, then multiplying both 
sides by 2, and adding 8 to both sides gives 


N N 

® + 2 IX - 2 E h m = 

m=l m-1 


N 

X_ V 2 

f 


AH 

F 

m=l 


or using (4.45) and reverting to the n subscript, 


N 


N 


(2N+l) 8- ~ ) a 2 = 1 - h Q - 2 ) h,. 


n 


(4.50) 


n=l 


n-1 


Multiplying both sides of (4.49) by m , summing from 1 to N, using 
(4.46) — (or 4.40 with j=2) — and reverting to the n subscript gives 

II 

*F - 


N N 

X X 1 4 X 1 2 v, 

n = - ) n h 

L. n 

n=l 


2f 

n=l s n=l 

Vie solve (4.50) and (4.5l) for 8 and X . 
Let 


H 


Qi . 1 - h o - 2 ; h n 


n-1 


N 


**2 ’ Z n \ 


n=l 


(4.51) 
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Then 


N 



n=l 


X = 


2f s [S lV (2N+l)Q 2 ] 

(2N+l)S 2 -2S 1 2 


and 


8 = 


Q l S 2 +2S l Q 2 

(2N+l)S 2 -2S 1 2 


Then 


h = h +8 
o o 


and from (4.49), for n > 1 


h = h + R 
n n ° 


2f 


= h + 
n 


Q l S 2 +2S l Q 2" n [S lV ( 2N+1)Q 2 ] 

( 21 ^ 1 ) 82 - 28 ^ 


Note that 


The constraint for p=l is obtained t>y letting X= 0 in (4 
Then we have 


h = 
n 


+8 


n 


n — 0,1, . . ., N 


«* 


(4.52) 


(4.53) 


(4.54a) 


(4.54b) 


• 54b). 
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where 


8 = 


l-h -2 
o 



n-1 

2N+1 


i-h n (° ) 

2N+1 * 


4.7 Band-pass filter 

The ideal band-pass smoothing filter transfer function is 


B I< f) 


= < 


f < f < f 
c — — c 


0 < f < f , f > f . 
— c 7 c 


( 4 . 55 ) 




f < 0 


Gee figure 4.7. 



Note that Bj(f) can be written as the difference of two ideal smoothing 
filter transfer functions [see ( 3 . 7 )] H^( f ) and II 2 (f), where H 2 (f) 
has cut-off f c and H^(f) has cut-off Then the weight function 

b (t) is 


b(t) = h 2 (t) - h^(t) 


(4.56) 
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where 


h x (t) < > H^f), 

h 2 (t) < > H 2 (f), 

and 

b(t) < > B^f). 

A useable design is obtained by taking the difference of two 
low-pass smoothing filters of the types discussed in section 4.2 and 
section 4.4. The difference of two Martin-Graham filters, each with 
roll-off length Af gives a satisfactory filter. The weight function 
of the resulting band-pass filter is then given by ( 4 . 56 ) with h^(t) 
and h 2 (t) the weights of the Martin-Graham filters. The weights of 
the corresponding digital filter are given by (3-44) and (4.56), 



* 

= J- th 2 (f£) -hj^)]. (4.57) 

s s s 

Now suppose B(f;f Q ) is a band-pass smoothing filter with the 
mid-points of the "pass bands" at + f Q , "pass band" width 2Af, and 
roll-off width Af. For purposes of illustration, we assume that 
B(f;f Q ) has the Martin-Graham type roll-off [see (4.18)]. Let 


0 < f < Af 


H(f) = i |[1+ co st f . -ffi . « 3 Af < f < Af+Af (4.58) 


0 


f > Af+Af 


(^H(-f) f < 0 . 

See figure 4.8. 
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Then H(f) is the transform of 


For f > 0 


and for f 


Thus 


1 XT- 

Af Af+Af 

Figure 4.8. 


_ sin 2itAf + sin 2«(Af+Af)t 
2«t( 1 - 4Af 2 t 2 ) 


B(f;f o ) = H(f - f o ), 

< 0 

B(f;f o ) = H(f + f Q ). 

B(f;f o ) = H(f - f o )+H(f + f Q ), 


see figure 4.9* 



Taking the inverse transform of each side and using the shift theorem 
(1.21) gives 


2«if t -2* if t 

b(t;t Q ) = h(t)(e ° + e ° ) 


= 2h(t)cos 2«f Q t. 


(4.61) 


The weights of the corresponding digital filter are given by 


f 

b (f ) = 2h cos 2*n~, 
no n r 


(*.&) 


where h = 4- h(^). 
nr r 
s s 

For a given f , the weights can be computed from (4.62) more 
o 

quickly than from (4.57)* If several successive filtering operations 
are to be performed for a set of f Q values, say f^f^, . . • ,f k , 
then, using (4.62), 

f, 

b n (fj) = 2h n cos 2njt^4 , j = 1,2, . . .,k. 

s 

But in order to use (4.57) the functions h^t) and h 2 (t) must be 

changed for each new value of f . and the entire expression must be 

J 

recomputed. 

From (4.62) we see that the error s' of a band-pass smoothing 
filter may be as much as twice the error e of the smoothing filter 
whose transfer function is H(f). 

In a manner analogous to the ideal smoothing case in section 3.2, 

the transfer function of a filter which will simultaneously "band-pass" 

"til 

filter and find the n — derivative is 


B n (f) = (2«if) n B(f) 


(4.63) 
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where B(f) is the transfer function of a band-pass smoothing filter. 
Then if b(t) « > B(f), 


b (n) (t) = b n (t) 


B n (f), 


and the weights [see the derivation of (4.3b)] are given by 

n A, 

b” = (- 1 ) (f. ) -4, 


where 




The weight function of a filter having several pass bands, each 
of equal pass width and roll-off width, can easily be found from 
(4.6l). Let + f , + f , . . . , + f be the mid points of the pass 

12 K 

bands, and denote the transfer function by B(f;f^,fg, . . • jf^) 

(see figure 


Figure 4.10. 


Then the weight function is 


b(t;f, ,f_, . . .,f ) = 2( V cos 2nf .t) h(t) 



The weights are given hy 



(4.67) 



CHAPTER V 


MARTIN-GRAHAM FILTERS 


5.1 Introduction 

We sh all call a filter a Martin-Graham filter if its transfer 
function either uses the Martin-Graham roll-off [see (4.l8)] or is 
derivable from a transfer function having the Martin-Graham roll-off. 

In section 4.2, we discussed the Martin-Graham smoothing filter 
and found its weight function h(t) [see (4.13)]* From h(t) and the 
formula (3.44) for computing the weights of the approximating 
digital filter, we found the weights h n [ see(4.15) ] which are used in 
the basic formula of digital filtering. 




( 3-45) 


where the g. are the input data values and the f . are the smoothed 
output values. A simple but detailed example of an application is 
given in Chapter VII. 

A Martin-Graham band-pass smoothing filter is easily obtained 
from the smoothing case and the discussion of section 4.7. 

In this chapter, we shall derive the weights for some other 
Martin-Graham filters. When referring to a set of data {®J ’ 
we assume that the data arises from a function g(t) such that 


1) g(t) = g(t) + p(t), where p(t) is a polynomial in t, 

2) g(t) satisfies conditions I - III of section 4.1, 

3) gjjj = g(]T-) where f g is greater than twice the highest 

s 

frequency in g(t). 


! 


76 


I 


Let g^ be the first data value and gjj be the last. If p(t) is not 

M M 

identically zero for ~ < t < then, in order to pass p(t) or 

s s 

differentiate it, constraints are necessary. Those for smoothing are 
in section 4.6. A general procedure is given in appendix A for the 
derivative cases, and the constraints for passing the first deriva- 
tive of p(t) will be given in the next section. 


5.2 Smoothing and first derivative filter 

We have shown that the transfer function of a filter which will 
smooth and find the first derivative to be 

r^w) = iwH(w) 

where H(w) is any smoothing filter transfer function [Put n=l and 
w = 2«f in (4.2)]. Note that Y^(w) inherits the cut-off, w , and 

V 

termination, w^, frequencies from H(w). 

Putting n=l in (4.3b), we obtain the weights of this filter in 
terms of the smoothing weights 


-f 


s 


dk 


(5.1) 


where h(t) < > H(w) 

The Martin-Graham 
terms of the frequency 



smoothing filter weights given by (4.15) in 
w f 

ratio, r = " = f~> are 


« 



cos knr, sin kit(2r +r ) 
d c a 

k«(l-4r d Sc^) 


sin 2nr T k + sin 2nr c k 
2rtk(l-4r d 2 k 2 ) 
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Then 


p p 

r T cos 2itr^k + r c cos 2itr c k 11^(1-12^ k ) 


1 * 
\ = - f s 


k(l-4r k ) 
a 


k(l-4r d "k c ) 


= -f 


r^ cos 2 Jtr,jk + r^ cos 2jtr c k -h^(l-12r d 2 k 2 ) 


k(l-4r d 2 k 2 ) 


(5-2) 


Note that y^ k = -y^, and by applying L' Hospital's rule, y^ = 0. 

In a manner analogous to that of section 4.6, we find that in 
order to pass exactly the derivative of P(t) of degree p the following 
conditions must be satisfied by the approximating filter transfer 
function 

N 

Y^(w) = 2i Y y^ sin |£. (5-3) 

n=l s 

(1) Tfj(0) = 0 


( 2 ) 


(3) 


Since 


dw 


= l 

w=0 


d P Y^(w) 

dw 5 

d^w) 

dw^ 


= 0 for p > 1. 


is odd for all even p > 0, (l) and (3) are automati- 


cally satisfied for even integers p > 0. In particular, if p(t) is 

of degree 2, we need to satisfy only (2). The constrained weights 

y are given by 
K 


- 1 1 

y k = y k + “^ > 


k > 1, 


(5.4) 


y k = - y k’ 
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where 


S. 




^ =Z ” 2 • 


(See [11] for the derivation, for p=4 from which the case p=2 follows 
easily. ) 

The constrained transfer function is 


V w) = 2i X y n sin F * (5*5) 

n=l S 

In order to smooth and differentiate a set of data u where the 

polynomial content is of degree 2 or less, put h n =y n in ( 3 . 45 ). 

This gives 


A K 

If we let h^ = 


and 



n=-N 


then 



V v) 






sin 


wn 

f 


s 


Then 



n=-N 


(5.6) 


(5-7) 


(5.8) 
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A 1 A 1 

Using 7 n = -y , we have 


N 


A1 


f m " f sl V n tg wn ‘Vn 1 ’ 
n=l 


(5.9) 


Writing 


in 

Y^w) = iwH(w) = wH(w)e 2 , 


we see that Y^w) has a phase shift of 90°. wH(w) is shown in figure 
5.1 for the Martin-Graham smoothing and first derivative filter. 



Other first derivative filters with different roll-offs have 

been examined and it was found that the Martin-Graham filter yielded 

the same or a more accurate result. 

In an attempt to avoid the lengthy computation of (5*l) for 

the weights y? - , a "three-point derivative" of the smoothing weights 
i£ 

h^ has been examined. Let 

8 - V 1 "Vl . (5-10) 

k 2_ 

f 

s 

With H(f) the transform of the weight function h(t) from which the 
are computed, we have 




[h k+l _h k-l ] 
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- | th(^i) -h (^ti) ] 


_s 

.2 


- I 1 » (f) 


-2*i(k+l)|- -2ni(k-l)|- 

[e s -e s ]df 


-f 

£ 

~2 

f 


■ I I H(f) 


2 -2jtik|- -2«i~ 2ni|- 


*[e -e s ]df 


-f 




2xif- -2«if- 


r e -e 


-f 


2i 


s “ 2rtik f 

-][-iH(f)e s ]df 


«k 




-f 


2 -sin 2«|- -2«ik|- 

-s [(2«if)H(f)]e s df. 

2n— 


The actual weights are 

f 

s 
2 


y£ = Y- f (2rtif)H(f) 


-2%±k— 

S -> ry 

e df . 


-f 


Comparing (5.11) and (5.12), we see that if we define a weight 
*y k = -0^, then the transfer function of the *y k is 


(5.11) 


(5-12) 


8i 



sin 2it- 


2rt|- 

s 


-(acif)H(f) 


which is the product of the desired transfer function and 


Now 


1 ~ 1 
y k - y k = F 


f 

sin 2«^r- 
s 
f 

2jij~ 

s 


_s f f 

2 sin 2n— -2rtik— 

f ( 2«if)H(f)[l jr-^le £ 


df 


2«: 


and 


_s 

2 


• o f 

sin 2n— 

1 i 0 


2rt- 


■f* 

for \-x-\ small. If the cut-off f is small, then H(f) in the above 
I c 

s 

integral becomes zero for relatively small. Then the are 

s 

good approximations of the y • It has been found empirically that 

f c 

for filters such that — < .1, the y^ give an acceptable output. 

s 

5-3 Band-pass smoothing and first derivative filter 

We have shown that the transfer function of a band-pass filter 
which will smooth and find the first derivative to be 


B (w) = iwB(w) 
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where B(w) is any hand-pass smoothing filter transfer function 
[put n=l and w=2xtf in (4. 63) ] • Note that B^"(w) has the same cutoff 
and termination frequencies as B(w). B(w) may be designed by either 
of the methods discussed in section 4.7. 

Putting n=l in (4.65), we obtain the weights of this filter 
in terms of the band-pass smoothing weights 


b 1 - f db * 
a = -f — — 

k s ^ 


where b(t) < > B(w) and \ b(^). 


s , s 


If the b k are obtained by taking the difference [see (4.57)] 
of the weights of two low pass filters, say h£ and h£, then 


b i = . f K <1 

k s 1 dk dkj * 


(5.13) 


When the b^ are obtained by the second method [see (4.62)], 
we have 


> 4 - 


-2f 


= -2f 


-( 


cos 2k«- 


( 5 . 14 ) 


-h 2jtf f dh 

k —sin 2k«^2 + -^cos 2k«- 
s s 


f dh f 

= 4nh, f sin 2k«^ -2f -#cos 2k*-£ . 

K o f s dk f 

s s 

To obtain a Martin-Graham filter of this type by the first method, 
we simply select two Martin-Graham smoothing filters with transfer- 

weight functions h'(t) < »H*(w) and h"(t) < > H"(w) and compute 

the weights b^ by (5.13)» To use the second method, the appropriate 
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Martin-Graham filter with h(t) < — > H(w) is selected and the weights 
are computed by (5.14) . These weights are used for the in 
(3.45). Note that a factor of f can be removed from the sum (3- ^5) 
in a manner analogous to the first derivative case [see ( 5«7) a > n{ i 


(5-8)]. 



(a) Smoothing transfer functions 





Y^w) - Y X (w) 


Figure 5*2. 


5.4 Smoothing and second derivative filter . 

Letting n=2 and w=2rtf in (4.2), we find that the transfer 
function of a filter which will smooth and find the second derivative 
is 

Y 2 ( w ) = -w^H(w ) (5.15) 
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where H(w) is any smoothing filter transfer function. 

Putting n=2 in (4.3b), we find that the weights of the filter 
in terms of the smoothing weights are 




(5.16) 


where h(t) <E H(w) and h^, = j- h(^). 


Using the Martin-Graham smoothing weights given by (4.15) in 

terms of r = and ( 5 . 1 6 ) gives 

s s 


f 2 


y k = 


k(l 


2y 

'^ r 2 k 2 ) [2kr i\- -^p(l-12 r ^ 2 )-2K(r 2 sin 2*^ 


+r^ sin 2«r c k)] ( 5 * 17 ) 


where y fc is given by equation (5.2). 

For k=0, using L' Hospital' s rule gives 


y 2 = f2 je(r T +r c ) . li|_( r 3 +r 3jj . (5-18) 

This gives the weights to be used in the formula (3.45). Note that 
a factor of f^ may be removed in this case. 

A constraint is developed in appendix A to improve the fit of 
the approximating transfer function at some specific frequency 
ratio r. 
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CHAPTER VI 


INTEGRATION 


6.1 The Romberg technique . 


One approach to the problem of integrating a set of numerical 
data is to smooth the input data and then apply one of the standard 
numerical integration techniques. A recently developed technique 
which seems to be well suited to the smoothed output of a filter is 
the one developed by Romberg (see appendix B). If 


I = 


a+2 k At 



ds 


( 6 . 1 ) 


k k 

where f(s) is known at the 2 +1 points: a,a+ j At, 0=1,2, ...,2, 

At = ■j— , then we approximate I by the diagonal element in the 
s 

4* Vi 

k— row of the Romberg array, i.e., 


I i 2 k AtT _. 

k,0 


( 6 . 2 ) 


6 . 2 Integrating filters . 

Another approach to the problem more in keeping with the one of 
this report is to design a filter which will simultaneously smooth 
and either give the indefinite integral or a definite integral of 
the input function. 

PlTT ff. 

Let Ae be a component of an input to a filter. Assuming 

that the constant of integration is zero, the indefinite integral 

of this component is (2rtif) If this is to be the output 

of the filter, then, using (3-5) > we find that the transfer function 
must be 

x(f) = (2«if)" 1 (6.3) 
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Letting k'(t) = g(t) and assuming that the constant of integration is 
zero and that k(t) satisfies conditions sufficient for the Fourier 
integral theorem to hold, we have 

k(t) <■ > (2rtif) _1 G(f) (6.4) 

where g(t) < > G(f). 

If we also smooth, we have 

k(t) < > (2nif) _1 H(f)G(f) ( 6 . 5 ) 

where H(f) is the smoothing filter transfer function. Then the 
transfer function of a filter which will simultaneously smooth and 
give the indefinite integral is 

y^“^(f) = (2rtif) _1 H(f). (6.6) 

Note that the smoothed output, g(t), of the smoothing filter is the 
inverse transform of H(f)G(f) and that 



g(a)dD£. 


(6.7) 


For the transfer functions, H.(f), 0=1,2, . ..,5, of the smoothing 

filters discussed in Chapter IV, Y _1 ^(f) has an infinite discontinuity 
at f=0. Hence, in order to approximate Y^”^(f) with a truncated 
Fourier series, we must modify Y^”^(f) on an interval containing zero. 
To avoid some integrals which cannot be evaluated in closed form, we 
shall consider only the case j=l, i.e., an Ormsby type filter. 

Let Af > 0 and 

|f(Af)" 2 , | f | < Af , 


Y (- 1 )(f) 


= (aci) 


-1 


,-l 


Af < f < f , 
— - c 5 


( 6 . 8 ) 


f -f 
T 

f Af 
c 


f c < f < f T > 

If I > f T , 


0 
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and Y^C-f) = -Y^ _1 \f) for f < 0. See figure 6.2. 



f • Af 

The weights in terms of the frequency ratio r = — , r^ = — > 
f f s s 

r c “ r» r T = r> are 

s s 


y ( - X) = 2 
n 2jt r 


f cos 2; 

2* 2 r,f l n 

d s v * 


2njtr , sin 2rutr 
d d 


2itr ,n 
a 


r cos 2n«r (sin 2n«r - sin 2rutr ) 


T 


r n 
c 


2«r n 
c 


- 2nr d [Si (2nnr c ) - Si(2rutr d )] 




where 


siW -/ HT* ^ = f (gillla-i) • 

0 k=l 


1) k.i i[ (a-i) 


Define Y <_1 ^ 0. Also, Y*" 1 ^- -Y^" 1 ' 
o * -k k 

A definite integral 

From ( 6 . 7 ) we have 


t+a 


k(t+a) - k(t-a) = J g(a)dD£, 


t-a 


( 6 . 9 ) 
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and by ( 6 . 5 ) and the shift theorem. 


k(t+a) - k(t-a) < > (2jtif) -1 G(f)H(f) [e 2rtiaf _ e " 2niaf ] 

or 

k(t+a) - k(t-a) < > H(f)G(f). (6.10) 

Thus, if H(f) is a smoothing filter transfer function, the transfer 
function of a filter which will simultaneously smooth and give the 
integral of the input over [t-a,t+a] is 




Let 


( 6 . 11 ) 


x(t) = « 


then 


t < a 


|t| > a , 


( 6 . 12 ) 


x(t) < > X(f) = 


sin 2naf 
nf 


Applying the convolution theorem gives 

co 

y^ _1 ^(t) = f h(z)x(t-z)dz 
-00 


(6.13) 


t+a 

= f h(2 


h(z)dz 
t-a 


(6.14) 


- f h(t. 

-a 


z)dz 


where y^ _1 ^(t) < > * X \f) and h(t) < > H(f), 
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CHAPTER VII 


APPLICATIONS 

7*1 Editing and det erminat ion of digital filter parameters . 

In order to apply a digital filter to a set of data K) ’ 
we assume that the data values are obtained by talcing equally spaced 
samples of a function g(t) which satisfies the three conditions of 
section 5.1. A variety of problems may arise from the methods used 
to obtain the samples, and editing may be necessary. Common problems 
are missing values and "bad" values, i.e., values grossly in error. 
Since these can affect the output considerably, it is important to 
replace them in some manner. The common practice is to consider the 
"bad" values as missing values and then replace each missing value 
by linear interpolation between the nearest data values on each 
side of the missing value. (See [7])* 

Next, the following parameters must be determined: 

A. The largest frequency, f , which Is present in the data. 

This is commonly found by visually determining the shortest 
period in the data. 

B. The sampling frequency, f , which must be at least 2f . 

S 

C. The cut-off frequency, f^, which is chosen to be at least 
as great as the highest frequency of interest present in 
the data. 

D. The termination frequency, f^. This should be chosen such 

that either, (l) no frequencies present in the data sire in 

the interval (f ,f ) or, (2) frequencies appearing in (f ,f ) 
0 T c X 

have no significant amplitude. 

E. The value of N and hence the number of weights, 2N+1, of the 
filter. 
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From the above, the corresponding frequency ratios may be found 

f 


from r = 


That is, r c = — , 
s 


fl 

r T f * 
s 


Af 

r d f 


7.2 Empirical error bounds ror Mart in-Graham filters 

Empirical error bounds are found by recovering the digital 
filter's transfer function, i.e., computing 

N 2nrti^A 

vv =Zv 

n=-N 


j=l,2, . . . ,k, for various values of the parameters of section 7*1* 

The recovered values are then compared with the designed or ideal 

transfer function values at the f.. An expression for the error € 

J 

is then determined in terms of N and the other parameters. 

The following error bounds were obtained by transfer function 
recoveries and comparison with bounds obtained by the method of 
section 4.5. 

I. Martin-Graham smoothing filter. 

For a maximum error e [see (3*47)1 of about .01, take 


N > 


1*25 
r , 


1.25f 
s 

Af 


(7*1) 


This gives a maximum error of 1% (e referred to unity) between 
the actual transfer function and the designed transfer function. 
Note that the error does not change with r c , r^ held constant. 

The bound given by the method of section 4.5 was compared 
with the results of computation with r c values ranging from 
.025 to .2, r^ values ranging from .021 to .11, and N values up 
to 100. It was found to be about 5 times too large. Hence, in 
terms of the frequency ratio, 


e 


1. 

5* ° e 4N 5_? 


r d - 1 


where "log" denotes the natural logarithm. 


(7*2) 
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II. Martin-Graham first derivative filter. 

Comparison of recoveries for f > f^, i.e., where Y^(f) 
is ideally zero, and the bound obtained by the method of section 
4.5 yielded, over the same range of frequencies ratio and N 
values given above, the expression 


€' 


f 

= -f Ur c +r T ) log 


4N 2 r 2 -l 
a 


«N(4N 2 r d 2 -l) 


(7-3) 


III. Martin-Graham second derivative filter. 

As above, the following expression was found 

f 


2 lafir 2 


2 2 
Wr +r„) log 


€ " “ ' 2 l "'*c"T / - w& ^2 2 


r +r 
c T 


r d - 1 


N(4N 2 r d 2 -l) 


]. 


(7-4) 


IV. Martin-Graham band-pass filters. 

The error can be as much as the sum of the errors in the 
low-pass filters from which the band-pass filter is derived (see 
section 4.7). Hence, in band-pass smoothing the error may be 
twice that obtained with a low-pass smoothing filter having the 
same roll-off length Af. 

The values of e' given by (7-3) become too large for small r d , 

but are still useable for r =.021. The values of e" given by (7*4) 

d 

are too small for large r, and small r . The actual value may be as 
4 d c 

much as — e" for r, values from .07 to .11 and r values of .025 to 
3d c 

.07. However, it is still useable, e' and e" are values for the error 

on the rejection band |f| > f ( |r| > r^) . The error on the pass -band 

|f| < f ( |r| < r ) is essentially the same. For the first derivative 

filter, the amplitude at f ideally is 2nf c = 2jtf s r c . For an error 

of 1 % of 2nf , we need 
c* 
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e' = .0l(2*f r ) 
s c 

= (.02)«r f 
c s 

= (. 08 )(«r c )-| . 

Comparing with ( 7 . 3 ), we see that N must he taken such that 


4nV 2 


( .08 ) k' c = (r « T ) log -5-3- ♦ ■ ( 7 .S) 

a & 

2 2 

For the second derivative, the amplitude at f ideally is 4n f Q = 

4* 2 f Similar to the above, we find that for an error of 1 % 

s c 

2 2 

of f , we need to take N such that 
c ’ 


. 2 2 
4N r . 
d 


( . 08 )* 2 r 2 = «(r 2 + r 2 ) log ~ ■ f - + 

c C T 4N 2 r ,-l N(4Hr. -l) 

d d 


(7-6) 


7.3 Programs and examples for Mart in-Graham filters 

When the appropriate filter has been chosen, the data edited, 
and the parameters of section 7*1 determined, the filtering can 
be performed. The weights, h^, of the filter are computed from the 
weight function h(t) and (3*44). If the data has a polynomial 
content, constrained weights h R are computed from the h n (see 
section 4.6, section 5.2, and appendix A). Then the output is 
computed using (3*45)* 

Consider the following examples. We take as the input the 
function 

g(t) = A-^cos 2-sif^t + AgSin 2«f 2 t 

+ A^cos 2rtf^t + A^. (7*7) 
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Using the Mart in-Graham filters, we perform the operations of smoothing; 
smoothing and finding the first derivative; and smoothing and finding 
the second derivative. The sampled version of g(t) is 

= g(f-) = AjCOS 2*2- f + A 2 sin 2*§-f 
s s s 

+ A cos 2*j-f + A, , (7.8) 

0 s J 

f 

apd going to the frequency ratio, r = — , gives 

s 

g n = A^os 2nnr 1 + A 2 sin 2rcnr 2 + A^cos 2rtnr^ + A^. (7*9) 

The following PDQ FORTRAN program was run on the IBM 1620. It is 
sectioned by comment cards which state what each part of the program 
does. Table 7*1 gives the frequencies used for the various runs. Table 
7-2 gives the frequency ratios, coefficients of the input and desired 
outputs, the corresponding program symbols, and the values of these 
parameters for each run. The value of N used, and hence the number of 
weights for each run, is given by the last two digits in the run num- 
ber. That is. Run 2.20 reads, "Run 2 with N = 20". 

The following symbolism was selected: 

FS: The sampling frequency, f . 

s 

HO: The central smoothing weight, h Q . 

DDHO: The central smoothing and second derivative weight . 

H(l): The smoothing weights h^, 1^0. 

DH(l): The smoothing and first derivative weights. 

DDH(l): The smoothing and second derivative weights, 

I 4 0 . 

TF1: The recovered transfer function for smoothing. 

TF2: The recovered transfer function for smoothing and 

the first derivative divided by 2n. 

TF3: The recovered transfer function for smoothing and 

2 

the second derivative divided by . 

Z(l): The input samples g R , n = I-N-l. 

Zl(l): The input samples on the range of interest. 

Z2(l): The desired smoothed output. 
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DZ(l): The desired smoothed and first derivative 


output . 


TlTY7 ( T \ « 

JUJULJ \-LJ 


The desired smoothed and second derivative 


output. 

Si: The actual smoothed output. 

S2: The actual smoothed and first derivative output. 

S3: The actual smoothed and second derivative output* 

The following weight properties are used in the program: 

1) Smoothing: h_ n = h n » 

1 1 

2) Smoothing and first derivative: y_ n = -y n 

2 _ 2 

3) Smoothing and second derivative: y_ n - Y n • 


4) Smoothing and integrating: y 


= v<:“= -v<-» 


The results for each run follow the program. 

In Run 2.20, the error e in the transfer function is about .011. 
At t = .8, the actual output differs from the desired output by .037* 
Referred to the desired output, this is an error of about 2%. 
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C C WEIGHTS. RECOVERY. FILTERING OF DATA. PDQ FORTRAN 
DIMENSION A ( 30 ) » H(30), DH<30). DDHI30). Z(lOl) 
DIMENSION ZlUO). Z2<40)» DDZ<40), DZ(40) 

P=3. 14159 
READ 1.XN1.TC.TD 

1 FORMAT (3F10.0) 

READ 5.T1.T2.T3.FS 

5 FORMAT (4F10.0) 

READ 14.A1.A2.A3.A4 

14 FORMAT (4F10.0) 

READ 18.B1.B2.B3.B4 
18 FORMAT (4F10.0) 

N1=XN1 

C C UNCONSTRAINED WEIGHTS 
H0=2.*TC+TD 
RT = T C+TD 

DDHO=TD*TD* ( RT+TC ) / ( P*P*P )-(RT*RT*RT+TC*TC*TC ) / ( 6.*P ) 
DDH0=8.*P*P*P*DDH0 
DO 2 1=1 ♦ N 1 
X = I 

A ( I >=SINF(HO*X*P)/ <X*P) 

H ( I ) =A( I )*COSF<TD*X*P)/< 1,-4.*Td*Td*X*X) 

DH ( I )=RT*C0SF<RT*X*2.*P)+TC*C0SF<TC*X*2.*P) 

DH( I ) =DH ( I ) -H( I )*< l.-12»*TD*TD*X*X ) 

DH( I ) = DH ( I ) / ( X* ( l.-4.*TD*TD*X*X) ) 

DDH ( I ) =-2.*DH( I )*( l.-12.*TD*TD*X*X )+24.*TD*TD*X*H( I ) 
DDH ( I ) =DDH ( I )-2.*P*TC*TC*SlNF(2.*TC*P*X) 

DDH ( I ) =DDH ( I )-2.*P*RT*RT*SlNF(2.*RT*P*X ) 

2 DDH ( I ) =DDH ( I ) / ( X* ( 1 . -4 ,*TD*TD*X*X ) ) 

N2 = N1 

C C CONSTRAINED SMOOTHING WEIGHTS 
S1 = 0 

DO 3 1=1. N2 

3 S 1=S 1+H ( I ) 

FN1=2*N1+1 

S 1 = 1 •-( HO+2 »*S1 ) 

S1=S1/FN1 
DO 4 1=1, N1 

4 H ( I ) = H ( I ) +S 1 
HO=HO+Sl 

C C TRANSFER FUNCTION RECOVERY 
Z1=100«*TC-1. 

Z2=Zl+6. 

Z 3 = Z2 + 100,*TD 
Z4=Z3+6. 

DO 13 K=1 .57 
TF1 = 0 
TF2=0 
TF3 = 0 
X =K 

IF(X-Zl) 6,66.66 

6 Y=K-1 
Y=.01*Y 
GO TO 10 


66 

IFCX-22) 

7,77,77 

7 

Y = Y + * 005 



GO TO 10 


77 

I F ( X-Z3 ) 

8,88,88 

8 

Y=Y+.01 



GO TO 10 
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88 I F ( X-Z4 ) 9.99.99 
9 Y=Y+.005 
GO TO 10 
99 Y = Y + . 0 1 

10 CONTINUE 

DO 11 1=1. N2 
X=I 

TF1=TF1+2.*H< I )*COSF< 2.*X*P*Y) 

TF2 = TF2 + 2.*DH(I)*5INF( 2.*X*P*Y ) 

11 TF3=TF3+2.*DDH< I ) *COSF ( 2 . *X*P*Y ) 

T F 1 = HO+T F 1 

TF2=-FS*TF2/ <2.*P) 

TF3=DDHO+TF3 

TF3=FS*FS*TF3/(4.*P*P) 

Y 1= Y*FS 

PUNCH 12.Y1.TF1.TF2.TF3 

12 FORMAT (E7.3* 3E15.8) 

13 CONTINUE 

C C INPUT DATA 
P2=P*P 
M1=N1+1 
M=2*Nl+40 
DO 15 1=1, M 
T = I -M 1 

C l = COSF ( 2.*P*T1*T ) 

S=SINF(2.*P*T2*T> 

C3=COSF(2.*P*T3*T ) 

15 Z (I) =A1*C1+A2*S+A3*C3+A4 

C C INPUT DATA ON THE RANGE OF INTEREST 
DO 17 1=1,40 
T = I - 1 

C1=C0SF(2.*P*T1*T) 

S=S INF ( 2 .*P*T2*T ) 

C3 = C0SF< 2.*P*T3*T ) 

Z 1 ( I ) =A1*C1+A2*S+A3*C3+A4 
Y=T/FS 

PUNCH 16 » Y »Z 1 ( I ) 

16 FORMAT (E7.3, E15.8) 

17 CONTINUE 

C C DESIRED OUTPUTS 
DO 20 1=1,40 
T = I - 1 

C 1 = COSF < 2»*P*T1*T ) 

S=SINF(2.#P*T2*T ) 

C3=COSF( 2.*P*T3*T ) 

Z 2 ( I )=B1*C1+B2*S+B3*C3+B4 

DDZ. ( I ) =-4.*P2*FS*FS* ( Bl*Tl*Tl*Cl+B2*T 2*T2*S+B3*T 3*T3*C3 ) 
C1 = SINF(2.*>P*T1*T) 

S = C0SF< 2.*P*T2*T ) 

C3=SINF<2.*P*T3»T) 

DZ ( I )=-2.*P*FS*(Bl*Tl*Cl-B2*T2*S+B3*T3*C3) 

Y = T/FS 

PUNCH 19,Y,Z2( I) ,DZ< I ) »DDZ( I ) 

19 FORMAT (E7.3, 3E15.8) 

20 CONTINUE 

C C ACTUAL OUTPUTS 
DO 23 K= 1 » 40 
M2=K-1 
M3=Nl+l 
S1 = 0 
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S2=0 
S3=0 
T =M2 
T=T/FS 

DO 21 I = 1 * N 1 
K1=M3-I 
K2= I+M2 
K3=M3+I+M2 

S1=S1+H( K1 ) *Z ( K2 )+H( I ) *Z ( K.3 ) 
S2=S2-DH ( K1 )*Z(K2)+DH< 1 )*Z(K3 ) 

21 53=S3+DDH(K1)*Z ( K2 ) +DDH ( I )*Z( K3) 
K4= M3+M2 

S l=HO*Z ( <4 ) +S1 

S2=-FS*S2 

S3=DDHO*Z(K4)+S3 

S3= FS#F S*S3 

PUNCH 22 »T *S1 »S2»S3 

22 FORMAT ( E7 . 3 » 3E15.8) 

23 CONTINUE 
STOP 
END 


TABLE 7.1 
FREQUENCIES 


FREQ. 

RUN 1.20 

RUN 2.20 

RUN 3.20 

RUN 4 < 

f l 

.5 

.5 

.5 

.5 

f 2 

.9 

.95 

o 

. 

csi 

2.0 

f 3 

2.0 

1.9 

4.0 

4.0 

f 

10.0 

10.0 

10.0 

10.0 

s 

f 

o 

• 

*— 1 

1.0 

2.0 

2.0 

c 

Af 

.6 

• 6 

• 6 

. 6 


1.6 

1.6 

2.6 

2.6 
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TABLE 7.2 

PROGRAM SYMBOLS AND PARAMETER VALUES 


FREQ. 

RATIO 

PROGRAM 

SYMBOLS 

RUN 1.20 

PARAMETER VALUES 
RUN 2.20 RUN 3.20 RUN 4. 

r 

c 

TC 

.1 

.1 

.2 

.2 

r T 

RT 

.16 

.16 

.26 

.26 

r d 

TD 

.06 

.06 

.06 

.06 

r l 

T 1 

.05 

.05 

.05 

.05 

r 2 

T 2 

.09 

.095 

.2 

.2 

r 3 

T 3 

.2 

.19 

.4 

.4 

COEFFICIENTS* 





A 1 

A 1 

1.0 

1.0 

1.0 

1.0 

A 2 

A2 

1.0 

2.0 

2.0 

2.0 

A 3 

A3 

.5 

1.5 

1.5 

1.5 

A 4 

A4 

.5 

1.0 

1.0 

1.0 


Bl 

1.0 

1.0 

1.0 

1.0 


B2 

1.0 

2.0 

2.0 

2.0 


B3 

0.0 

0.0 

0.0 

0.0 


B4 

.5 

1.0 

1.0 

1.0 


*The coefficients Bl, B2, B3, and B4 are for the desired output. 
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Run 1.20. Recovered transfer functions 


f TF1 TF2 

.000E-50 .10000000E 01 -.OOOOOOOOE-50 

. 100E 00 .99797082E 00 .10644201E 00 

. 200E 00 .99514616E 00 .20274402E 00 

• 300E 00 .99689773E 00 .29448908E 00 

.400E 00 .10028117E 01 .39483205E 00 

,500E 00 .10056242E 01 .50387711E 00 

, 600E 00 .10002252E 01 .60704498E 00 

, 700E 00 .99285677E 00 .69780276F 00 

,750E 00 .99266489E 00 .74289948E 00 

.800E 00 .99624218E 00 .79162209E 00 

, 850E 00 .10028210E 01 .84539667E 00 

.900E 00 .10095384E 01 .90275637E 00 

.950E 00 .10114761E 01 .95858893E 00 

.100E 01 , 1 002 3041 E 01 • 10044468E 01 

.110E 01 .92563225E 00 .10252742E 01 

. 1 2 OE 01 . 75029065E 00 .90216597E 00 

. 1 30E 01 . 50077482E 00 .64488536E 00 

• 140E 01 .25002918E 00 .34476281E 00 

.150E 01 .72924460E-01 .11260134E 00 

( 160E 01 -.32426400E-02 . 19521540E-02 

, 165E 01 11377890E-01 -. 13822814E-01 

, 170E 01 84137 400 E-02 -. 1.3544579E-01 

.175E 01 11392600 E-02 -. 57 59352 1 E-02 

. 1 80E 01 . 5260 3800E-02 .26600947E-02 

• 185E 01 . 79699400E— 02 . 76060798E-02 

.190E 01 .65757700E-02 . 78909192 E-02 

,200E 01 -,21774900 E-02 373921 19 E-04 

« 210E 01 -.58878800E-02 52 53 1985 E-02 

• 220E 01 -» 58409 000E-03 -. 16588389E-02 

, 2 30E 01 .46315200E-02 , 32 1 3030 5 E-02 

, 240E 01 .26071700E-02 . 24427648 E-02 

, 250E 01 -.27621400E-02 14786985E-02 

, 260E 01 -.36313500E-02 24903026E-02 

,270E 01 .64339000E-03 , 12064591 E-03 

, 280E 01 .35876200E-02 . 20251410E-02 

, 290E 01 , 12612300E-02 , 77 1 31372 E-03 

,300E 01 -,26257000 E-02 1 3044897E-02 

,310E 01 -.25406600E-02 -, 11997365E-02 
» 320E 01 , 10930 100E-02 , 54057674E-03 

,330E 01 .29542400E-02 , 12258077E-02 

. 340E 01 ,55387 000E-03 , 96642464E-04 

, 350E 01 -»24929200E-02 -, 96748445 E-03 
.360E 01 -.18729000E-02 -.49995495E-03 
, 370E 01 , 13618800E-02 .58565360E-03 

,380E 01 ,25360l00E-02 , 6647430 1 E-03 

• 390 E 01 .76130000E-04 -, 2C482956E-03 

,400E 01 -.24060800E-02 -,63557844E-03 
,410E 01 -, 13945 100E-02 -.90762957E-04 
, 420E 01 .15722800E-02 .48410550E-03 

, 430E 01 .22277100E-02 .25444567E-03 

,440e 01 -, 30863000 E-03 -.30625590E-03 
,450E 01 -,23597500 E-02 -.3 1051454E-03 
.460E 01 -.10061 100 E-02 . 1 594269 1 E-03 

, 470E 01 , 1772 1900E-02 .30706919E-03 

,480E 01 , 19841 000E-02 -.58141577E-04 

, 490E 01 -.65 48 5000 E-03 -.28422375E-03 


. 500E 01 -.23 43 3100 E-02 -. 19021419E-07 


TF3 

-. 131187 59E-02 
-. 10946376E-01 
-. 39 3517 67E-01 
-.87443287E-01 
-. 1 5869925E 00 
-.25 324265E 00 
-. 3641 1447E 00 
-, 48710996E 00 
-, 55597912E 00 
-.63277556E 00 
-.71 8 9 7448 E 00 
-.81 362836E 00 
-.91 209207E 00 
-, 10055591E 01 
-, 11271631E 01 
-, 1 081 1058E 01 
-.83823027E 00 
-.48412609E 00 
-, 1 6974264E 00 
-.20212799E-02 
, 24411391E-01 
, 24458559E-01 
. 10723431E-01 
-.523138 1 IE-02 
-. 15449837E-01 
-. 16755633E-01 
-, 40204820 F.-03 
. 12614394E-01 
. 49615192E-02 
-.837321 44E-02 
-.77437910E-02 
. 37839675E-02 
,855 1 17 53E-02 
, 7 0627570E-03 
-.731 7 3091E-02 
-.42105583 E-02 
. 46011 1 16E-C2 
.61 619848E-02 
-.11957794E-02 
-.64062478 E~02 
-, 2 1 880954E-02 
, 49 185261E-02 
. 46417541E-02 
-.22487438 E-02 
-. 5 59347 07 E-02 
-.74151273 E-03 
, 50422546E-02 
, 3403 5565E-02 
-.31038074 E-02 
-.49991779 E-02 
,3361 1 586E-03 
, 50543043E-02 
.23483579E-02 
-.37513775 E-02 
-.43462630 E-02 
. 14058413E-02 
.51270327E-02 


101 



Rim 1.20. Input over the range of interest 


t 


Zi(t) 


.000E- 

-50 

.20000000E 

01 

. 100E 

00 

•21413919E 

01 

. 200E 

00 

• 18093361E 

01 

.300E 

00 

• 16753913E 

01 

« 400E 

00 

• 17 340 388E 

01 

.500E 

00 

.13090206E 

01 

.600E 

00 

.96808960E- 

■01 

«700E 

00 

-•12212563E 

01 

# 800E 

00 

-• 16958 1 32 E 

01 

. 900E 

00 

-•1 2 263299E 

01 

« 1 OOE 

01 

- • 5877 89 10E 

00 

# 1 1 OE 

01 

-• 35933910E 

00 

• 120E 

01 

-.23177496E 

00 

.130E 

01 

.38400309E 

00 

.140E 

01 

.13435080E 

01 

.150E 

01 

• 18090 172E 

01 

• 160E 

01 

•13316611E 

01 

« 170E 

01 

•49590509E 

00 

« 180E 

01 

•2199589 4E 

00 

» 190E 

01 

•63697291E 

00 

• 200E 

01 

• 1 04894Q6E 

01 

* 2 1 OE 

01 

• 968 1 4 547E 

00 

• 220E 

01 

• 7 79 1 7 473E 

00 

• 2 30E 

01 

• 1 1 090 443E 

01 

• 240E 

01 

• 1 8078 412E 

01 

.250E 

01 

• 2 0000 066E 

01 

.260E 

01 

•11898459E 

01 

.270E 

01 

- #66488 190E- 

-01 

.280E 

01 

-.83885020E 

00 

.290E 

01 

-.93397380E 

00 

.300E 

01 

-.95105200E 

00 

« 3 1 OE 

01 

-.12651218E 

01 

«320E 

01 

- • 1 398 0 790E 

01 

« 330E 

01 

-•67970750E 

00 

• 340E 

01 

• 71357 536E 

00 

• 350E 

01 

• 18089978E 

01 

.360E 

01 

•19615601E 

01 

.370E 

01 

• 1 5595957E 

01 

• 380E 

01 

•13862599E 

01 

• 390E 

01 

• 1 542 7700E 

01 
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Run 1.20* Desired outputs 


t 

• 000E-50 

Smoothed 

Z2(t) 

.15000000E 

01 

1— Derivative 
DZ(t) 

•56548620E 01 

2— Derivative 
DDZ(t) 

-.98695878E 01 

• 100E 

00 

• 19868829E 

01 

•38037554E 

01 

265209C6E 

02 

• 200E 

00 

.22 138439E 

01 

• 561 14911E 

00 

-•36918729E 

02 

• 300E 

00 

.20799007E 

01 

— •32503322E 

01 

-.37526523E 

02 

.400E 

00 

• 15795 324E 

01 

-•65923646E 

01 

-.27688979E 

02 

.500E 

00 

.80902060E 

00 

-.85196792E 

01 

-.988 16658E 

01 

.600E 

00 

-.57702570E-01 

-•84650394E 

01 

•11002238E 

02 

« 700E 

00 

-.81675000E 

00 

-•641 26361E 

01 

•29111678E 

02 

• 800E 

00 

-.12913022E 

01 

-.29062224E 

01 

• 39395685E 

02 

• 900E 

00 

13808338E 

01 

• 11108 59 IE 

01 

• 39 1 18474E 

02 

• 1 OOE 

01 

10877891E 

01 

•45748552E 

01 

•28665593E 

02 

• HOE 

01 

-.51385310E 

00 

•66144971E 

01 

•11394594E 

02 

• 1 20E 

01 

• 17272979E 

oo 

• 680 19809E 

01 

-• 7 42041 76E 

01 

. 130E 

01 

•78851564E 

00 

•52658747E 

01 

-•22 2 207 44E 

02 

. 140E 

01 

.11890066E 

01 

• 26327922E 

01 

-• 2 8 864472 E 

02 

• 150E 

01 

•13090172E 

01 

-•18222171E 

00 

-• 2 5870408E 

02 

• 160E 

01 

• 1 177 1446E 

01 

-• 22699080E 

01 

-• 14821747E 

02 

• 170E 

01 

•90040828E 

00 

-.30130995E 

01 

•19056161E 

00 

• 180E 

01 

• 62447 306E 

00 

-.22756563E 

01 

. 13905253E 

02 

. 190E 

01 

•48247401E 

00 

-.43553584E 

00 

•21586242E 

02 

• 200E 

01 

. 54894060E 

00 

• 17474137E 

01 

• 205428 83E 

02 

• 210E 

01 

.81362636E 

00 

•33863211E 

01 

• 10996902E 

02 

• 220E 

01 

• 1 1836763E 

01 

• 37636986E 

01 

-•39765130E 

01 

• 230E 

01 

• 1 5 135 599E 

01 

•25751092 E 

01 

-. 19416283E 

02 

• 240E 

01 

• 16533449E 

01 

•42256899E- 

-01 

-• 30049203E 

02 

« 250E 

01 

• 15000066E 

01 

-•31415225E 

01 

-• 31977530E 

02 

• 260E 

01 

.10353242E 

01 

-.60178027E 

01 

-• 23949877E 

02 

• 270E 

01 

•33801 188E 

00 

-.76582527E 

01 

-.78 145864E 

01 

• 280E 

01 

-.43433296E 

00 

-♦ 74568799E 

01 

• 11992048E 

02 

• 290E 

01 

-• 10884676E 

01 

-•53280230E 

01 

• 29769379E 

02 

• 300E 

01 

-• 14510 520E 

01 

-• 17475506E 

01 

• 4028 18 18E 

02 

• 310E 

01 

-• 14196459E 

01 

•23770048E 

01 

•40359512E 

02 

« 320E 

01 

-•99358050E 

00 

• 5968 7 149 E 

01 

•29875160E 

02 

• 330E 

01 

-•27518875E 

00 

•80962657E 

01 

• 11793730E 

02 

• 340E 

01 

.55908404E 

00 

• 82456119E 

01 

-. 87212573E 

01 

• 350E 

01 

• 13089978E 

01 

•64655110E 

01 

-• 2 5869906E 

02 

• 360E 

01 

• 18070335E 

01 

• 33430089 E 

01 

-. 349641 09E 

02 

.370E 

01 

•19640926P 

01 

-• 18254309E 

00 

-.33823464E 

02 

« 380E 

01 

• 17907802E 

01 

-• 3 1087390E 

01 

-• 23390364E 

02 

• 390E 

01 

• 13882812E 

01 

-.467 28742 E 

01 

-•73792113E 

01 



Run 1.20* Actual outputs 


Smoothed 1 — Derivative 2 — Derivative 


t 


SI 


S2 


S3 


. GOOE- 

■50 

• 15045 354E 

01 

.56721806E 

01 

-*1003 1 434E 

02 

• 100E 

00 

• 1 997 0064E 

01 

.38109615E 

01 

-* 2 67477 57E 

02 

. 200E 

00 

• 22279052E 

01 

.55427308E 

00 

-.37171376E 

02 

• 300E 

00 

.20935504E 

01 

-# 327 22 799 E 

01 

-.37763351E 

02 

.400E 

00 

• 15882 834E 

01 

-.66266840E 

01 

-*278672 56E 

02 

. 500E 

00 

• 8 1 08 7 941 E 

00 

-.85605104E 

01 

-.99597611E 

01 

. 600E 

00 

-.62 149062E- 

■01 

-.85048700E 

01 

.11049072E 

02 

« 700E 

00 

-.82612818E 

00 

-.64441304E 

01 

* 29 27 18 69E 

02 

.800E 

00 

-.13043409E 

01 

-.29238554E 

01 

*39620 48 2E 

02 

« 900E 

00 

-.13953876E 

01 

• 11095947E 

01 

♦39345057E 

02 

. 100E 

01 

-•11001086E 

01 

.45888662 E 

01 

* 2 8 8439 74E 

02 

« 1 10E 

01 

-, 5201 3754E 

00 

.66394216E 

01 

* 1 1 496990E 

02 

• 120E 

01 

* 1 7 36 5 562 E 

00 

.6831 5450E 

01 

-.74053338E 

01 

. 130E 

01 

*7944491 8 E 

00 

• 52938569E 

01 

-.22290500E 

02 

.140E 

01 

• 1 1964 5 1 8 E 

01 

.26547608E 

01 

-* 2 89962 24E 

02 

• 150E 

01 

* 1 3 1 56452 E 

01 

-. 16804116E 

00 

-.26020127E 

02 

« 160E 

01 

♦11820574E 

01 

-.22627312E 

01 

-.14942387F 

02 

# 170E 

01 

* 90 28 0 760 E 

00 

-•30103346E 

01 

. 1 2 26 8 360E 

00 

. 180E 

01 

*62 33747 5E 

00 

-. 22740332 E 

01 

. 13880260F 

02 

.190E 

01 

*4782 47 76E 

00 

-.43242794E 

00 

. 21574888E 

02 

.200E 

01 

* 54440 436E 

00 

•17527655E 

01 

•20517271E 

02 

. 2 10E 

01 

•81255887E 

00 

. 3392 2 50 1 E 

01 

. 109381 09 E 

02 

• 220E 

01 

•11879120E 

01 

. 37666320E 

01 

-.408161 1 6 E 

01 

• 2 30E 

01 

.15218 07 3E 

01 

.25710047E 

01 

-* 19571985E 

02 

• 240E 

01 

• 1662 8 000E 

01 

.28258741E- 

-01 

-* 30238055E 

02 

• 2 50E 

01 

• 1 5084 564E 

01 

-•31658828E 

01 

-•32 1546 16E 

02 

• 260E 

01 

* 1041 3034E 

01 

— .60 50 140 5 E 

01 

-.24059611E 

02 

. 270E 

01 

• 33964794E 

00 

-. 76935647E 

01 

-.78 197993E 

01 

• 280E 

01 

-.43919722E 

00 

-. 748845 14E 

01 

. 12 09408 IE 

02 

• 290E 

01 

-.11002328E 

01 

-•53490081E 

01 

*299540 85E 

02 

.300E 

01 

-♦1466 8 36 5 E 

01 

-.17529024E 

01 

*405122 43E 

02 

« 310E 

01 

-• 14345700E 

01 

.23889526E 

01 

. 40591663E 

02 

• 320E 

01 

-.10037791E 

01 

• 59957294E 

01 

* 3 0057299E 

02 

• 330E 

01 

-•27940 1 5 1 E 

00 

•81329167E 

01 

. 11876349E 

02 

.340E 

01 

.56052073E 

00 

•82847717E 

01 

-.87627748E 

01 

. 3 50E 

01 

•13156258E 

01 

.65000505E 

01 

-* 26019624E 

02 

• 360E 

01 

* 18 179545E 

01 

. 33673742 E 

01 

-* 35174977E 

02 

• 370E 

01 

• 19766 374E 

01 

-• 171 1 122 IE 

00 

-* 34043700E 

02 

. 380E 

01 

•18008065E 

01 

-.31096643E 

01 

-* 2 3 5 82 42 OE 

02 

« 390E 

01 

* 13926950E 

01 

-.46827420E 

01 

-.752031 55E 

01 


See figures 7*1 and 7.2, pages 114 and 115 


Run 2.20. Input on the range of interest 


t 


Zi(t) 


• 000E- 

-50 

« 35000000E 

01 

• 100E 

00 

•36274107E 

01 

• 200E 

00 

•25751 186E 

01 

• 300E 

00 

•21823776E 

01 

• 400E 

00 

.2 7722949E 

01 

• 500E 

00 

•27394576E 

01 

• 600E 

00 

. 79557440E 

00 

• 700E 

00 

-•20318853E 

01 

• 800E 

00 

-•32932433E 

01 

• 900E 

00 

-• 1904419QE 

01 

• 100E 

01 

•59547310E 

00 

• 1 1 OE 

01 

• 18734 1 52E 

01 

. 120E 

01 

•14509463E 

01 

• 1 30E 

01 

•92990740E 

00 

• 140E 

01 

• 16398419E 

01 

• 150E 

01 

•27896500E 

01 

. 160E 

01 

•25112433E 

01 

• 170E 

01 

•45319630E 

00 

. 180E 

01 

-• 14425950E 

01 

• 190E 

01 

-• 10865011E 

01 

• 200E 

01 

• 12879099E 

01 

• 210E 

01 

•33852534E 

01 

« 220E 

01 

•35193545E 

01 

• 2 30E 

01 

• 23964943E 

01 

• 240E 

01 

•18789242E 

01 

• 250E 

01 

•24141995E 

01 

• 260E 

01 

•24604273E 

01 

.270E 

01 

•64480020E 

00 

• 280E 

01 

-•21 362846E 

01 

• 290E 

01 

-.34471 109E 

01 

• 300E 

01 

-•20816204E 

01 

« 3 10E 

01 

•52717530E 

00 

• 320E 

01 

•20028098E 

01 

• 3 30E 

01 

• 1 7244565E 

01 

• 340E 

01 

•12223375E 

01 

• 350E 

01 

•19002987E 

01 

• 360E 

01 

•30762402E 

01 

• 370E 

01 

•28730387E 

01 

• 380E 

01 

•81532020E 

00 

• 390E 

01 

-• 12359836E 

01 
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The transfer functions for Run 2.20 are the same as for Run 1.20. 
The frequencies f^ and f^ were chosen such that they appear where the 
transfer function has a relatively large error for each in smoothing. 
Run 2.20. Desired outputs. 


t 


Z2(t) 


DZ(t) 


DDZ(t) 


,OOOE- 

■50 

• 2 000000QE 

01 

•11938042 E 

02 

-# 98695878E 

01 

. 100E 

00 

• 30752 22 5E 

01 

• 8 9 0 2 92 24E 

01 

-.49439679E 

02 

« 200e 

00 

•36685695E 

01 

• 2 548 1 1 86E 

01 

-.74239047E 

02 

• 3 00E 

00 

.3539620 IE 

01 

-•51437837E 

01 

-.75343517E 

02 

.400E 

00 

.2678 1 152E 

01 

-.1 169 02 70 F. 

02 

-.31829733E 

02 

• 500E 

00 

.13128752E 

01 

-• 14932650E 

02 

11 147463E 

02 

« 600E 

00 

16056850E 

00 

-•137897 10E 

02 

• 33 390021E 

02 

.700E 

00 

- . 13092 641 E 

01 

-.86185972E 

01 

.67136173E 

02 

« 800E 

00 

18050697E 

01 

-• 10970380E 

01 

,79 102480E 

02 

,900E 

00 

-• 15 313713E 

01 

.63460558E 

01 

.65691925E 

02 

.100E 

01 

-.61804350E 

00 

•11353726E 

02 

. 3 1 8 89993E 

02 

. 110E 

01 

.60691420E 

00 

. 12434841E 

02 

104935^1E 

02 

. 12 OE 

01 

•17319999E 

01 

•94562217F 

01 

-.96920600E 

02 

. 130E 

01 

• 24033 347E 

01 

• 366 5 140 2 E 

01 

-.65140906E 

02 

.140E 

01 

• 2443 5 997E 

01 

-.27632956E 

01 

59394571E 

02 

.150E 

01 

• 19079905E 

01 

- • 749 5 242 5 E 

01 

-# 32351089E 

02 

.160E 

01 

. 10583625E 

•01 

-•88560850E 

01 

• 58806474-E 

01 

.170E 

01 

•26517030E 

00 

-.64133147E 

01 

• ^-1 322436E 

02 

» 180E 

01 

128 14780E 

00 

-.1 12 23806E 

01 

.61034916E 

02 

• 1 90E 

01 

.69286900E-01 

•5014580 7E 

01 

.37659391E 

02 

,200E 

01 

• 82441 320E 

00 

• 96 5802 43E 

01 

. 32015644E 

02 

« 2 1 OE 

01 

•18882153E 

01 

• 1096 1359 E 

02 

-.7 1 475097E 

01 

• 220E 

01 

♦ 28806 549E 

01 

•82331297E 

01 

-.46 166189E 

02 

,230E 

01 

• 34232 900E 

01 

• 2 1997 10 3E 

01 

-.71 198657E 

02 

« 240E 

01 

. 32736022E 

01 

-• 5 224642 7E 

01 

-.73046333E 

02 

.250E 

01 

• 24 142 373E 

01 

-• 11582958E 

02 

5 0 387988E 

02 

• 260E 

01 

• 10657772F 

01 

-.14714394E 

02 

10303582E 

02 

.270E 

01 

-.38205070E 

00 

-*• 13497869 E 

02 

. 34100402E 

02 

.280E 

01 

-*• 1497 6 5 3 5E 

01 

- • 8 24346 1 8 E 

01 

.68 149563E 

02 

• 290E 

01 

19500681 E 

01 

-• 5960 18 1 8E 

00 

• 80609807E 

02 

• 300E 

01 

16180 5 18E 

01 

• 70 1 68343 E 

01 

,675 1 9496E 

02 

« 3 1 OE 

01 

-.62 856450E 

00 

• 1 2 2 02990 E 

02 

. 33525547E 

02 

, 320E 

01 

.68832620E 

00 

. 13409593E 

02 

-. 97354072E 

01 

.330E 

01 

.19124074E 

01 

. 104 36504E 

02 

47649667E 

02 

• 340E 

01 

.267 51995E 

01 

• 4484262 3E 

01 

-.67646423E 

02 

.350E 

01 

•27820194E 

01 

-• 2277985 IE 

01 

634921 83E 

02 

.360E 

01 

• 2272 5463E 

01 

-• 747 34449 E 

01 

37379897E 

02 

.370E 

01 

• 13995972E 

01 

- • 9 3434638 E 

01 

. 90358681E 

00 

• 380E 

01 

• 534192 30E 

00 

- • 7 3 5 1 9577 E 

01 

• 37436183E 

02 

• 390E 

01 

• 30476800E- 

-01 

-.2359996 2E 

01 

. 59042124E 

02 
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Run 2 , 20 * Actual outputs 


Smoothed 1 — Derivative 2 — Derivative 


t 


SI 


S2 


S3 


, 000E- 

-50 

•20154878E 

01 

•12045973E 

02 

11041620E 

02 

. 100E 

00 

.30971031E 

01 

.89155165E 

01 

-• 5040419 IE 

02 

• 2 00E 

00 

• 3687 2 694E 

01 

•25226223E 

01 

-.74375279E 

02 

• 300E 

00 

•36564000E 

01 

-• 5 1 57371 IE 

01 

-.75311883E 

02 

» 400E 

00 

.2696 1 842E 

01 

-.11717893E 

02 

-.52501832E 

02 

« 500E 

00 

.13258467E 

01 

-• 1 504063QE 

02 

12261396E 

02 

«600E 

00 

-.16579159E 

00 

-• 1 3967839E 

02 

•33067783E 

02 

« 700E 

00 

-.13370774E 

01 

-.87584182E 

01 

« 68 2 89526E 

02 

« 800E 

00 

18423 126E 

01 

10952590E 

01 

•80894529E 

02 

.900E 

00 

-.15573091E 

01 

• 64767 1 32 E 

01 

• 66607086E 

02 

» 100E 

01 

-.62278052E 

00 

. 11500090E 

02 

•31397538E 

02 

• 110E 

01 

.61629707E 

00 

. 1 2 50 6 166 E 

02 

11472647E 

02 

« 120E 

01 

.17432870E 

01 

•94662856E 

01 

-#4726647 5E 

02 

. 130E 

01 

.24131902E 

01 

• 368 1 0689 E 

01 

-.64896809E 

02 

. 140E 

01 

• 24566893E 

01 

-.27293315E 

01 

-.59538814E 

02 

. 150E 

01 

.19242 0 84E 

01 

-• 75068831 E 

01 

-# 33 32992 IE 

02 

• 160E 

01 

.10667777E 

01 

-•89584918E 

01 

•49231610E 

01 

. 170E 

01 

.25453440E 

00 

-• 65483543E 

01 

•41571858E 

02 

. 180E 

01 

15447 227E 

00 

-• 11707356E 

01 

• 62 48261 OE 

02 

. 190E 

01 

.45440250E- 

•01 

• 5 1 0607 1 3 E 

01 

• 58962975E 

02 

• 200E 

01 

.81959411E 

00 

• 98 16Q733E 

01 

. 31974414E 

02 

• 210E 

01 

.19026873E 

01 

• 1 1 0 66 38 1 E 

02 

-.82875370E 

01 

• 220E 

01 

.29017031E 

01 

• 82426470 E 

01 

-.47149863E 

02 

• 2 30E 

01 

.34409083E 

01 

• 2 1686487E 

01 

-• 7 1 341 566E 

02 

• 240E 

01 

•32887141E 

01 

-• 52406639E 

01 

-.72959075E 

02 

• 250E 

01 

.24304671E 

01 

1 16Q9264E 

02 

-•50975316E 

02 

« 260E 

01 

•10775113E 

01 

-•14816202E 

02 

-.11380296E 

02 

• 270e 

01 

-.38771926E 

00 

-•1 36 70845 E 

02 

• 3 3745373E 

02 

• 280E 

01 

-.15257821E 

01 

-.83829090E 

01 

* 69 2632 33E 

02 

• 290E 

01 

-. 19882024E 

01 

-* 59 549004E 

00 

.82427016E 

02 

• 300E 

01 

-. 16452934E 

01 

•71310027E 

01 

•68515082E 

02 

• 3 10E 

01 

-.63408873E 

00 

• 123 5947 5 E 

02 

•33087545E 

C2 

• 320E 

01 

•69812756E 

00 

• 1349262 5 E 

02 

-•10741499E 

02 

• 3 30E 

01 

.19250822E 

01 

• 104538QOE 

02 

-• 4 8 0699 99E 

02 

• 340E 

01 

, 26866783E 

01 

# 450 2455 7E 

01 

-•67448987E 

02 

« 3 50E 

01 

.27966715E 

01 

-# 22424527E 

01 

-.63635556E 

02 

.360E 

01 

•22906269E 

01 

-♦ 7482060 IE 

01 

-.38367769E 

02 

• 370E 

01 

• 14104324E 

01 

-•9445 1399 E 

01 

-•12686000E 

00 

.380E 

01 

.52596078E 

00 

-• 74938559E 

01 

•37577600E 

02 

• 390E 

01 

.54566600E- 

-02 

-• 242 24347E 

01 

•60433601E 

02 



Run 3.20. Recovered transfer functions 


f 


.000E-50 

. 1 00e 

00 

.200E 

00 

• 300E 

00 

. 400E 

00 

, 500E 

00 

. 600E 

00 

, 700E 

00 

.800E 

00 

. 900E 

00 

, 100E 

01 

. 110E 

01 

• 1 20E 

01 

. 1 30E 

01 

. 140E 

01 

. 150E 

01 

. 160E 

01 

. 170E 

01 

. 175E 

01 

. 180E 

01 

.185E 

01 

, 190E 

01 

. 195E 

01 

• 200E 

01 

. 2 10E 

01 

, 220E 

01 

.230E 

01 

.240E 

01 

. 250E 

01 

.260E 

01 

. 265 E 

01 

. 270E 

01 

. 275 E 

01 

. 280E 

01 

. 285E 

01 

. 290E 

01 

. 300E 

01 

.310E 

01 

,320e 

01 

. 330E 

01 

. 340E 

01 

. 3 50E 

01 

. 360E 

01 

, 370E 

01 

.380E 

01 

, 390E 

01 

• 400E 

01 

.410E 

01 

. 420E 

01 

.430E 

01 

• 440E 

01 

.450E 

01 

.460E 

01 

.470E 

01 

. 480E 

01 

.490E 

01 

.500E 

01 


TF1 

.10000000E 01 
.99945238E 00 
.99863960E 00 
.99897643E 00 
* 10005Q81E 01 
• 10015378E 01 
• 10006253E 01 
.99870627E 00 
.99833411E 00 
• 1 0003624E 01 

• 10022 588E 01 
.10010699E 01 
.99788668E 00 
.99730807E 00 

• 10010904E 01 

• 1 0042 013E 01 

• 1 000792 7E 01 
.99441370E 00 
.99386400E 00 
.99666469E 00 
.10023570E 01 
• 10084103E 01 
• 1 0 101403E 01 
.10012704E 01 
.92604993E 00 
.75147569E 00 
.50108251E 00 
•24911331E 00 
,72 100960E-0 1 

-.28700700E-02 
-•10526310 E-01 
-.74207000E-02 
38653000 E-03 
, 54924300 E-02 
.76056200E-02 
•57678000E-02 
-•28776700 E- 02 
-• 55 256200 E-02 
•31000000E-03 
,48 150 100E-02 
• 18 449 600 E-02 
-.33 996600 E-02 
-, 32646400 E-02 
• 14945900E-02 
•37360000E-02 
•49900000E-03 
-•32386200E-02 
-•21445 3 00 E-02 
• 19569700E-02 
•30798600E-02 
-.25244000E-03 
-,31156 400 E-02 
-•14207000E-02 
•22879500E-02 
•26338200E-02 
-.83840000E-03 
-.30766900E-02 


TF2 

-,00000000e-50 
• 10761383E 00 
• 20400454E 00 

• 29440003 E 00 
•39296784E 00 
• 50204665E 00 
,608 35788E 00 
•70226832E 00 
•79240672E 00 
•89362575E 00 
• 10048168E 01 
•11093631E 01 

• 1 1994607E 01 
•12893453E 01 
, 13955016E 01 
•15101747E 01 
, 16095504E 01 

• 169 12238E 01 
, 17347289E 01 

• 17859937E 01 
• 18460814E 01 
• 19107547E 01 

• 19697331E 01 
•20077460E 01 
• 19529443E 01 
• 16525654E 01 

• 1 1436934E 01 
. 5936 1745E 00 
, 18685445E 00 
.61257277E-03 

-.24190123E-01 
-,22251017 E-0 1 
-.83783146E-02 
, 5 8 50691 2 E-02 
, 13754990E-01 
. 13636174E-01 
-,57320315 E-03 
-.94010325E-02 
-.29603702E-02 
,56742181 E-02 
, 458460 17E-02 
-. 2282 1846E-02 
-.45695297E-02 
-.39263350E-03 
,33800911 E“02 
. 19979871E-02 
-. 16358196E-02 
-,24325171 E-02 
-.24102604E-04 
. 19017655E-02 

• 1 1 19 1 1 18 E-02 
-.81849770E-03 
-, 14057018E-02 
-.30348199E-03 

.92689450E-03 
, 99433503 E-03 
.48763747E-07 


TF3 

.86044627E-03 
-, 10389035E-01 
-.41397524E-01 
-.89287873F-01 
-. 1 5734707E 00 
-, 25000005E 00 
-.36397262E 00 
-,49 198182E 00 
-. 63 5447 48E 00 
-, 804881 89E 00 
-. 10034423E 01 
-.12 186314E 01 
-, 14398548E 01 
-, 16784461E 01 
-, 19546465E 01 
-, 22632507E 01 
-.25729319E 01 
-, 28758868E 01 
-.30381955E 01 
-, 32 179358E 01 
-, 34178864E 01 
-. 36314467E 01 
-, 38398498E 01 
-.40 125424E 01 
-, 40982205E 01 
-. 36371321E 01 
-, 26347093E 01 
-. 14255807E 01 
-.46301103E 00 
, 18476202E-02 
.64760126E-01 
, 5 75 2692 1E-Q 1 
, 18153645E-01 
-. 2 1730061E-01 
-.42844089E-01 
-.39896803 E-0 1 
.74476768 E-02 
.32943817E-01 
. 54 166142E-02 
-.25297890E-01 
-. 14953360E-01 
. 1 563 16 55E-0 1 
, 20224603E-01 
-.48185396E-02 
-.207311 59E-01 
-.5494 2011 E-02 
. 1661 1028E-01 
. 13 303063E-0 1 
-.91247478 E-02 
-. 17237979E-01 
.99852194E-04 
. 1675 1 383E-0 1 
, 84910334E-02 
-. 12012786E-01 
-. 14552608E-01 
.43299174E-02 
, 16693073E-01 



i 


Run 3 * 20 . Input over the range of interest 


t 


Zi(t) 


« 000E-50 

• 3 5000000E 

01 

• 100E 

00 

•26396454E 

01 

« 200e 

00 

•344811Q7E 

01 

• 300E 

00 

.87375520E 

00 

.400E 

00 

- « 18 066 307E 

01 

« 500E 

00 

• 24999 908E 

01 

»600E 

00 

• 1 3795 795E 

01 

.700E 

00 

.20513030E 

01 

« 800E 

00 

-.52102270E 

00 

.900E 

00 

-•30667170E 

01 

• 100E 

01 

• 14999 788E 

01 

♦ 110E 

01 

• 73754350E 

00 

. 120E 

01 

•18300614E 

01 

. 130E 

01 

-.29977140E 

00 

. 140E 

01 

-.24246944E 

01 

« 1 50E 

01 

• 24999 642E 

01 

. 160E 

01 

• 19976 199E 

01 

. 170E 

01 

•32268552E 

01 

# 1 80E 

01 

• 10970 545E 

01 

« 190E 

01 

-• 1 1646317E 

01 

• 200E 

01 

•34999576E 

01 

• 2 10E 

01 

•26396713E 

01 

• 220E 

01 

•34480875E 

01 

• 230E 

01 

•87585440E 

00 

• 240E 

01 

-• 18066762E 

01 

• 2 50E 

01 

• 24999 5 36E 

01 

• 260E 

01 

• 13796 088E 

01 

• 270E 

01 

•20512811E 

01 

.280E 

01 

-•52092470E 

00 

• 290E 

01 

-• 30667660E 

01 

• 300E 

01 

• 14999 363E 

01 

• 3 10E 

01 

•73756630E 

00 

# 320E 

01 

• 1 8 300 321 E 

01 

.330E 

01 

-•29968060E 

00 

• 340E 

01 

- • 24247 30 1 E 

01 

• 350E 

01 

• 24999 164E 

01 

.360E 

01 

• 1 9976 393E 

01 

• 370E 

01 

• 32268247E 

01 

• 380E 

01 

•10971464E 

01 

• 390E 

01 

-.11646838E 

01 
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Run 3»20» Desired outputs 




Smoothed 


st ... 

1 — Derivative 

2— Derivative 

t 


* \ 

L,d\T,) 


DZ(t) 


DDZ(t) 


,000E- 

50 

, 2 OOOOOOOE 

01 

.25132720E 

02 

-.98695878E 

01 

,100E 

00 

.38531690E 

01 

.67956594E 

01 

-.30975557E 

03 

200E 

00 

.29845913E 

01 

-.22 179345E 

02 

-. 19362355E 

03 

, 3 00E 

00 

.41222060E 

00 

22874442E 

02 

. 17983632E 

03 

,400E 

00 

-.59309770E 

00 

.47785073E 

01 

• 29731968E 

03 

i 5 00E 

00 

, 99999080E 

00 

.21991130E 

02 

.16654732E- 

02 

,600E 

00 

.25930937E 

01 

.47787586E 

01 

-. 29731867E 

03 

,700E 

00 

.15877988E 

01 

-.228 7 4291 E 

02 

-. 17983906E 

03 

,800E 

00 

-.98457240E 

00 

-.22179509E 

02 

. 19362082E 

03 

, 900E 

00 

-.18531747E 

01 

.6795 3974E 

01 

. 30975660E 

03 

100E 

01 

-.21200000E- 

•04 

.25132712E 

02 

. 98729430E 

01 

110E 

01 

• 1951048 /l E 

01 

.87375126E 

01 

-« 29098146E 

03 

.120E 

01 

« 1366572. JE 

01 

-. 18486037E 

02 

-. 17765692E 

03 

.130E 

01 

-.76333620E 

00 

-. 17791409E 

02 

,19 143604E 

03 

140E 

01 

12111427E 

01 

. 10753908E 

02 

• 30342050E 

03 

. 1 50E 

01 

.99996420E 

00 

.28274310E 

02 

. 5 07 16949E- 

-02 

160E 

01 

.32111154E 

01 

. 107 54677E 

02 

-. 30341731E 

03 

170E 

01 

.27633813E 

01 

-. 17790923E 

02 

-. 19 144412E 

03 

180E 

01 

. 63 347450E 

00 

-. 18486488E 

02 

. 17764882E 

03 

. 190E 

01 

.48929400E- 

-01 

. 87367750E 

01 

. 29098460E 

03 

« 200E 

01 

.19999576E 

01 

.25132736E 

02 

986287 88E 

01 

« 2 1 OE 

01 

.38531574E 

01 

.67961828E 

01 

-. 30975352E 

03 

• 220E 

01 

.29846287E 

01 

-, 221790 17E 

02 

-. 19362901E 

03 

. 230E 

01 

.41225920E 

00 

-.22874747E 

02 

. 17983085E 

03 

• 240E 

01 

-.59310580E 

00 

.47780053E 

01 

. 29732171E 

03 

. 250E 

01 

.99995360E 

00 

.21991130E 

02 

• 83 22 5003E- 

-02 

.260E 

01 

• 25930 856E 

01 

.47792613E 

01 

-.29731664E 

03 

# 270E 

01 

♦15878375E 

01 

-• 22873987E 

02 

-. 17984453E 

03 

.280E 

01 

-.98453500E 

00 

-.22 179836E 

02 

. 19361536E 

03 

« 290E 

01 

18531862E 

01 

• 67948740E 

01 

. 30975866E 

03 

.300E 

01 

-.63700000E- 

-04 

.25132695E 

02 

. 98796524E 

01 

.310E 

01 

* 19510337E 

01 

• 87380040E 

01 

29097937E 

03 

.320E 

01 

.13666036E 

01 

-, 18485737E 

02 

-. 177662 31E 

03 

« 330E 

01 

-.76330610E 

00 

-.17791732E 

02 

. 19143065E 

03 

» 340E 

01 

-.1211 1609E 

01 

. 10753395E 

02 

.30342262E 

03 

. 3 50E 

01 

.99991640E 

00 

. 282743 10E 

02 

, 1 1833271E 

-01 

• 360E 

01 

.32110973E 

01 

. 10755190E 

02 

-. 303^15 19E 

03 

« 370E 

01 

.27634114E 

01 

-. 17790600E 

02 

-. 19 144950E 

03 

. 380E 

01 

.63350580E 

00 

-. 18486788E 

02 

. 17764342E 

03 

. 390E 

01 

.48914700E 

-01 

.87362830E 

01 

• 29098669E 

03 




Run 3»20. Actual outputs. 


Smoothed 


1— Derivative 


2— Derivative 


t SI S2 S3 

• 000E-50 .19966798E 01 .25230053E 02 -.88519600E 01 

•100E 00 .38609774E 01 .68308277E 01 -.31145924E 03 

• 200E 00 .29858278E 01 -.22280317E 02 -.19386771E 03 

,300E 00 .41012973E 00 -.22948935E 02 .18075636E 03 

• 400E 00 -.59110874E 00 .47872926E 01 .29749969E 03 

• 500E 00 .99513290E 00 .22075609E 02 .10193050E 01 

.600E 00 .25989650E 01 .48056716E 01 -.29902235E 03 

• 700E 00 • 15868873E 01 -.22978105E 02 -. 18008323E 03 

• 8 00E 00 -.98881099E 00 -.2225U55E 02 .19454088E 03 

.900E 00 -.18531235E 01 .68124424E 0l .30993660E 03 

•100E 01 64169700 E- 02 .25230049E 02 . 10890582E 02 

• HOE 01 . 19559322E 01 .87806289E 01 -.29268509E 03 

• 120E 01 • 13653205E 01 -.18571891E 02 -.17790110E 03 

• 1 30E 01 -.76723489E 00 -.17845091E 02 . 19235609E 03 

•140E 01 -.12101041E 01 . 10787154E 02 .30360045E 03 

•150E 01 .99510620E 00 .28384502E 02 , 10227100E 01 

• 160E 01 .32179365E 01 .10806049E 02 -.30512093E 03 

• 170E 01 ,27642 776E 01 -.17873931E 02 -. 19168832E 03 

• 180E 01 ,63172 394E 00 -.18343018E 02 . 17856888E 03 

,190E 01 .51905730E-01 .87617643E 01 .29116458E 03 

• 2 00E 01 , 19966371E 01 .25230066E 02 -.88452300E 01 

• 2 10E 01 .38609656E 01 .68313535E 01 -.31143716E 03 

.220E 01 .29858655E 01 -.22279989E 02 -.19387323E 03 

• 2 30E 01 .41016821E 00 -.22949241E 02 .18075090E 03 

.240E 01 -,59111 675E 00 .47867862E 01 .29750169E 03 

• 250E 01 .99509570E 00 .22075606E 02 .10259830E 01 

,260E 01 .25989569E 01 .48061764E 01 -.29902032E 03 

,270E 01 , 1586926 IE 01 -.22977801E 02 -. 18008878E 03 

,280E 01 -,98877 416E 00 -.22251484E 02 . 19453546E 03 

,290E 01 -.18 53 135 IE 01 .68119199E 01 . 30993869E 03 

.300E 01 -, 64594000 E- 02 .25230035E 02 .10897302E 02 

• 3 1 0E 01 , 19559174E 01 .87811226E 01 -.29268303E 03 

.320E 01 , 13653 520E 01 -.18571596E 02 -. 17790656E 03 

,330E 01 -.76720516E 00 -.17845418E 02 .19235080E 03 

.340E 01 -.12101223E 01 .10786646E 02 .30360260E 03 

• 350E 01 .99505880E 00 .28384506E 02 . 10294110E 01 

.360E 01 .32179183E 01 .10806554E 02 -.30511888E 03 

.370E 01 .27643073E 01 -.17873613E 02 -.19169365E 03 

,380E 01 .63175457E 00 -.18543311E 02 .17856364E 03 

.390E 01 .51891 750E-0 1 .87612840E 01 .29116652E 03 
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Run 4 . 


f 

.000E-50 
• 100E 00 
« 200E 00 
. 300E 00 
,400E 00 
• 500E 00 
,600E 00 
,700E 00 
. 800E 00 
.900E 00 

• 1 00E 01 
• 110E 01 
• 120E 01 
• 130E 01 
, 140E 01 
• 150E 01 
• 160E 01 
. 170E 01 
« 175E 01 

• 180E 01 
. 185E 01 

• 190E 01 
. 195E 01 
.200E 01 

• 2 1 OE 01 
. 220E 01 
.230E 01 
.2^0E 01 
• 230E 01 
.260E 01 
• 265E 01 
.270E 01 
. 275E 01 
. 28 OE 01 
» 285E 01 
.290E 01 
« 300E 01 
.310E 01 
.320E 01 
.330E 01 
.340E 01 
« 350E 01 
,360e 01 
« 3 70E 01 
.380E 01 
.390E 01 
• 400E 01 
.410E 01 
.420E 01 
.430E 01 
.440E 01 
« 450 E 01 
.460E 01 
.470E 01 
• 480E 01 
.490E 01 
, 5 OOE 01 


. Recovered transfer functions. 


TF1 

( 10000000E 01 
.99907998E 00 
.99931646E 00 
• 10010677E 01 
. 10000894E 01 
.99886238E 00 
.10006425E 01 
.10007752E 01 
.99880356E 00 
,99998 585E 00 
.10013162E 01 
.99909529E 00 
.99915694E 00 
, 10016514E 01 
.99977512E 00 
, 99814679E 00 
. 10017739E 01 

• 10010042E 01 
.99813180E 00 
•99640572E 00 
, 99803983E 00 
, 10026755E 01 

• 10060471E 01 
, 10005279E 01 
.92967818E 00 
.75327729E 00 

• 50035794E 00 
.24648207E 00 
.69765660E-01 
, 88390000 E-04 

•» 53850800E-02 
-» 25282200E-02 
, 14668500 E-02 
, 28774900 E- 02 

• 15368000E-02 
-.66490000E-03 
-» 12712400E-02 

, 1 1748000E-02 

• 17509000E-03 
-,9239 2000 E-03 

.40365000E-03 
.45518000E-03 
-.56294000E-03 
-.29990000 E-04 
.44878 00 OE-03 
-.22763 000 E-03 
-.22481000E-03 
.29683000E-03 
.18 100 000 E-04 
-.23712000 E-03 
. 10304000E-03 
. 13462000E-03 
-. 13365000E-03 
-.51090000E-04 
• 1 1976000E-03 
.97500000E-05 
-. 10946000E-03 


TF2 


00000000 E~ 

50 

• 10019297E 

00 

. 19971041 E 

00 

.30012350E 

00 

#4GG31356E 

00 

• 49948591 E 

00 

.60006845E 

00 

• 70066934E 

00 

• 79929887E 

00 

• 89969049 E 

00 

• 100 12 1 1 OE 

01 

• 10993790E 

01 

• 1 198866 3 E 

01 

.13018188E 

01 

# 1 4000287 E 

01 

• 1497 3674E 

01 

• 160238 13 E 

01 

* 1 70 1 8492 E 

01 

. 17471953E 

01 

#17 9 3940 4E 

01 

• 1 8463800E 

01 

• 190466 1 OE 

01 

# 196 12778 E 

01 

. 20008798 E 

01 

# 19528897E 

01 

• 16 570165E 

01 

. 11503241E 

01 

• 59209553 E 

00 

# 1743992 1 E 

00 

-.47688527E- 

-03 

- • 1 46 56341 E" 

-01 

~# 65433446 E* 

-02 

#47777932E“ 

-02 

• 86491 1 78 E- 

-02 

#43091 42 5E* 

-02 

-.26293922E- 

-02 

- • 39808054E- 

“02 

#45 2 1 8398 E- 

-02 

. 13739491E- 

-03 

-# 37339274E* 

-02 

•22937934E* 

-02 

• 1 7 0 2801 6 E' 

-02 

-.311 89940 E 

-02 

#52071151 E 

-03 

# 24772608 E 

-02 

-.2 1 068544 E 

-02 

89 52 526 5 E 

-03 

# 2 56959 1 3 E 

-02 

-.86277327E 

-03 

- # 1864442 6 E 

-02 

• 2077530 OE 

-02 

# 40070746 E 

-03 

-#22914146E 

-02 

• 1 1 5 34379 E 

-02 

. 14658183E 

-02 

-#21411032E 

-02 

18661849E 

-06 


TF3 

-.58318545 E-02 
-.80659397E-02 
-• 35478026 E-01 
-. 95055768 E- 01 
-. 16119351E 00 
-, 243958 72E 00 
36286209E 00 
-.49431150E 00 
-.63393498E 00 
-, 809568 36E 00 
-.10068540E 01 
-. 12056404E 01 
-. 14354893E Ol 
-. 169828 35E 01 
-. 19592377E 01 
-, 22406451E 01 
-• 2 5683347E 01 
-. 23955224E 01 
-. 30542592E 01 
-. 32228329E 01 
-. 34119223E 01 
-. 36206707E 01 
-.38304844E 01 
-, 40068302E 01 
-, 40957027E 01 
-, 36439244E 01 
-, 26522678E 01 
-, 14181303E 01 
-.43 186968E 00 
-.496783 66 E-02 
.32960773E-01 
. 17141369E-01 
-.77990086E-02 
-. 17510939E-01 
-.99314431E-02 
.35498696E-02 
.78814335E-02 
-.69 6 5 767 IE— 02 
-. 10527795E-02 
.55 163144E-02 
-.242882 48 E-02 
-.26208541E-02 
. 33 1840 10E-0 2 
.27863372E-04 
-.2 6C30469E-02 
. 13943591E-02 
. 11685898E-02 
-.17309233E-02 
.52383140E-04 
. 13160377E-02 
-. 64 149082E-0 3 
-.60460985E-03 
.76680001E-03 
. 16333002E-03 
-, 395997 53E-03 
-.33967984E-04 
.44758708C-03 
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The input and desired output for Run 4.30 are the same as for Run 
3.20. The number of weights was increased. Note the improved accuracy 
in the transfer function and the actual outputs. 

Run 4.30. Actual outputs. 


Smoothed 1 — Derivative 2— Derivative 


t 


SI 


S2 


S3 


.000E- 

-50 

• 1998 5248E 

01 

• 25 143768E 

02 

-•9 792 12 00E 

01 

.100E 

00 

*38 533630E 

01 

.68050311E 

01 

-.31032783E 

03 

• 200E 

00 

.29841874E 

01 

-* 22 1 94413E 

02 

-* 19395645E 

03 

• 3 OOE 

00 

.4 1 08 2 7 16E 

00 

-* 22 872747E 

02 

* 18008466E 

03 

• 400E 

00 

-.59418055E 

00 

* 47 800346E 

01 

* 29762004E 

03 

. 5 OOE 

00 

.99965360E 

00 

* 22005412E 

02 

1 59 366 00E 

00 

• 600E 

00 

•25947223E 

01 

* 47 902065 E 

01 

29819147E 

03 

« 700E 

00 

•15689840E 

01 

-• 2 2 888646 E 

02 

■“*18050508 E 

03 

.800E 

00 

-.98437667E 

00 

-.22178530E 

02 

. 1 93 5 3604E 

03 

.900E 

00 

-.18528241E 

01 

•67948506E 

01 

. 3 0975645E 

03 

. 100E 

01 

• 77930000E- 

-03 

• 2 5 1 43765 E 

02 

* 94734070E 

01 

.110E 

01 

•19534076E 

01 

* 8 7 448899 E 

01 

-* 29 200738E 

03 

• 120E 

01 

* 1368 0 091 E 

01 

-• 1 8 504907 E 

02 

-.17837570E 

03 

. 130E 

01 

- • 76339220E 

00 

-• 17794939E 

02 

* 1 9 140 3 99E 

03 

. 140E 

01 

- * 1 2 1 1 5 2 25E 

01 

* 10749293 E 

02 

•30357342E 

03 

. 150E 

01 

*999627 10E 

00 

*28282128 E 

02 

-* 1 5 596600E 

00 

. 160E 

01 

•32120402E 

01 

* 107 59978 E 

02 

30414267E 

03 

. 170E 

01 

*27632 290E 

01 

-.17810504E 

02 

-•19182978E 

03 

.180E 

01 

*63182935E 

00 

- * 1 848930 4E 

02 

* 17794993E 

03 

• 190E 

01 

•47116300E-01 

.87342318E 

01 

* 29 143807E 

03 

• 200E 

01 

* 19984823E 

01 

* 2 5 1 43779 E 

02 

-*9785 39 OOE 

01 

• 2 10E 

01 

.38533517E 

01 

.68055583E 

01 

-•31032575E 

03 

.220E 

01 

* 29842 2 50E 

01 

-•22194087E 

02 

-* 19396 191E 

03 

• 230E 

01 

* 4 1 086 573E 

00 

-• 22873052 E 

02 

. 18007914E 

03 

.240E 

01 

-• 59418 86 5E 

00 

•47795285E 

01 

* 297622 08 E 

03 

. 250E 

01 

•99961640E 

00 

* 2200541 1 E 

02 

15269400E 

00 

• 260E 

01 

• 25947 143E 

01 

.47907116E 

01 

- * 2 98 18948 E 

03 

.270E 

01 

* 1 5890228E 

01 

-.2 2 88 8340 E 

02 

-* 18051059E 

03 

. 280E 

01 

-.98433947E 

00 

-.221 78857E 

02 

•19353060E 

03 

.290E 

01 

-.18528357E 

01 

* 6794330 5 E 

01 

.30975857E 

03 

• 300E 

01 

* 73690 000E- 

-03 

* 2 5 14375 1 E 

02 

* 9480 1 2 00E 

01 

# 3 1 OE 

01 

• 19 533 929E 

01 

* 87 45 382 1 E 

01 

-* 29 200 533E 

03 

« 320E 

01 

* 1 3 68 0 404E 

01 

-.18504611E 

02 

17838113E 

03 

.330E 

01 

-• 76336249E 

00 

-.17795265E 

02 

* 19 139871E 

03 

• 340E 

01 

-.12115406E 

01 

• 10748786E 

02 

* 30 35 75 59 E 

03 

.350E 

01 

.99957970E 

00 

*282821 33 E 

02 

149272 00E 

00 

• 360E 

01 

•32120224 E 

01 

* 1 0760484E 

02 

-* 30414064E 

03 

• 370E 

01 

•27632586E 

01 

-.17810186E 

02 

19183507E 

03 

«38GE 

01 

•63 186026E 

00' 

-• 1 848959 5E 

02 

* 17794463E 

03 

« 390E 

01 

*47102350 E- 

-01 

* 8 7 3 37504E 

01 

* 29 144Q01E 

03 
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output 



7.4 Program and example for indefinite integration . 

The following program was run with the input given by (7*7) with 

A ± = 1.5, A ? = 2.0, A 0 = 1.5, Aj^ = o, f ± = 0.7, ? 2 = 0.9, = 2*0, and 

f =10.0. N = XN1 = 25. Also, 25 terms were used in computing the 
s 

sine integral — which is too many terms for small values of the argument. 

For large values of the argument, the first terms of the series may 

become large enough to cause loss of significance, and computation of 

the sine integral should be approached with caution. 

The other parameters chosen were f =1, Af = .6, Af the inner 

c 

and outer roll-off length. In terms of the frequency ratio, the input 
frequencies are .07, .09, and .2, and r Q = .1, r^ = .06, r T = .16. 

The notation chosen in the program is analogous to that of 
section 7.3. The results follow the program. In this case the 
weights are given along with the transfer function recovery, the 
desired output, and the actual output. 
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c c indefinite integration with smoothing 

DIMENSION TERM1(50),TERM2(50) ,A(30),b(30) 

DIMENSION C<30 J *H< 30) »Z< 101) *ZI (40 ) 

P = 3.14159 

READ 1. XM1 * XN1 * TC» TO 

1 FORMAT (FlO.CtFlO.OtFlO.OtFlO.O) 

READ 2. Tit T2 » T3t FS 

2 FORMAT (F10.0,F10.0tF10.0,F10.0) 

READ 3 1 Alt A2t A3 

3 FORMAT (FlO.OtFlO.OtFlO.O) 

READ 4, Bit B2t 83 

4 FORMAT (FlO.OtFlO.OtFlO.O) 

Ml = XM1 

N1 = XN1 
RT = TC + TD 
C C SINE INTEGRAL 
TERMl ( 1 ) =1. 

TERM2 ( 1 ) =1 • 

DO 8 I = It N1 

X = I 

XI = 2.*X*P*TC 
X2 = 2.*X*P*TD 
DO 5 K=1 tMl 

Y = K 
J = K + 1 

TERMl(J)=-Xl*Xl*TERMl(K)*(2.*Y-l.)/< < 2 .*Y+ 1 . ) * < 2 .*Y+1 . ) *2 .*Y ) 

5 TERM2 ( J ) =-X2*X2*TERM2(K)*(2.*Y-l,)/( ( 2.*Y+1. )*< 2.*Y+1. )*2.*Y ) 

51 = 0 

52 = 0 

DO 6 J = It Ml 
S1=S1+TERM1 ( J) 

6 S2=S2+TERM2 ( J ) 

C C WEIGHTS 

A ( I )=2.*P*TD*(X2*S2-X1*S1) 

B ( I ) = ( COSF (2.*P*X*TD)/X) - S INF ( 2 . *P*X*T D ) / ( 2 . *P*TD*X*X ) 
C(I) = (SINF(2.*P*X*RT) - SINF(2.*P*X*TC) ) / ( 2 . *P*TC*X*X ) 

H ( I ) = A ( I ) +B ( I ) +C ( I )-TD*COSF( 2.*P*X*TC) /( TC*X ) 

PUNCH 7t Xt H< I ) 

7 FORMAT ( F7 • 3 1 E20.8) 

8 CONTINUE 

C C TRANSFER FUNCTION RECOVERY 
DO 11 K = It 51 
HI = 0 

Y = K-l 

Y = Y * .01 

DO 9 I = It N1 
X = I 

9 HI = HI + 2.*H( I )#SINF( 2.*X*P*Y ) 

HI = H1/(2,*P*P*TD*FS) 

Y1=Y*FS 

PUNCH lO,YltHl 

10 FORMAT ( F7 . 3 1 E20.8) 

11 CONTINUE 

C C INPUT DATA 
M3=N 1+1 
M = 2*Nl+40 
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P 2 = P*P 

DO 14 I = 1 » M 

T = I -M3 

Cl = COSF(2.*P*Tl*T) 

S = SINF(2.*P*T2*T ) 

C3=COSF ( 2.*P*T3*T ) 

14 Z ( I )=A1*C1+A2*S+A3*C3 
C C DESIRED OUTPUT 

DO 16 I = 1» 40 
T = 1-1 

Cl = SINF(2.*P*T1*T) 

S = COSF ( 2 . #P*T2*T ) 

C3=SINF(2.*P*T3*T) 

Z I ( I ) = ( l./( 2.*P*FS) )*( (B1*C1/T1)-(B2*S/T2 ) + (B3*C3/T3) ) 
Y=T/FS 

PUNCH 15 *Y » Z I ( I ) 

15 FORMAT ( F7 . 3 1 E20.8) 

16 CONTINUE 

C C ACTUAL OUTPUT 

DO 19 K = 1* 40 
M4 = K-l 
M5 = Nl+1 
SI = 0 
T = M4 
T = T/FS 

DO 17 I = 1* N1 
K 1 = M5-I 
K2 = I + M4 
K3 = M5+M4+I 

17 SI = SI - H(K1)*Z(K2)+H( I )*Z( K3 > 

SI = S1/(2.*P2*TD*FS) 

PUNCH 18 » Tt. SI 

18 FORMAT ( F7 • 3 » E20.8) 

19 CONTINUE 
STOP 
END 
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Indefinite 

integration with smoothing. 

Weights and Transfer function. 

n 

Weights 

f 

Transfer function 



.000 

# 00 000000 E -99 

1.000 

-.21454526E-00 

.100 

-•48316595E-01 

2.000 

-. 3 5 92 465 5 E -00 

.200 

-.87799236E-01 

3.000 

-.39618853E-00 

.300 

-.12781742E-00 

4.000 

-•33380622E-00 

.400 

-•1790Q819E-00 

5.000 

-.21699964E-00 

.500 

-•22789864E-00 

6.000 

-.10052025E-00 

.600 

-.24757827E-00 

7.000 

-.22045076E-01 

.700 

-•23 157831E-00 

8.000 

. 10480938E-01 

.800 

-.20080206E-00 

9.000 

• 1 3045094E-0 1 

.900 

-•17572434E-00 

10.000 

• 8 5941740E-02 

1.000 

-•15644913E-00 

11.000 

• 1 1 353359E-01 

1.100 

-•13391 49 2E-00 

12.000 

•21432966E-01 

1.200 

-.10680382E-00 

13.000 

•30236526E-01 

1.300 

-.79328906E-01 

14.000 

•30126059E-01 

1.400 

-• 52227487E-0 1 

15.000 

•20557104E-01 

1.500 

-.25658076E-01 

16.000 

•73245800E-02 

1.600 

-.34776971E-02 

17.000 

-.28660550E-02 

1.700 

•22907104E-02 

18.000 

-.72928270E-02 

1.800 

•69471017E-03 

19.000 

-•78980180E-02 

1.900 

-•14145 40 7E-02 

20.000 

-.81692720E-02 

2.000 

-.27708135E-03 

21.000 

-. 9947546 OE-02 

2.100 

•99977962E-03 

22.000 

-.11789276E-01 

2.200 

• 16 964473E-03 

23.000 

-. 1 1421037E-01 

2.300 

-•76567173E-03 

24.000 

-.62803660E-02 

2.400 

-•13772748E-03 

25.000 

-.95407800E-03 

2.500 

• 60982346E-03 



2.600 

• 1 3 1 98459 E-03 



2.700 

-.49169692E-03 



2.800 

-•13306796E-03 



2.900 

• 39 5 33037E-03 



3.000 

• 1 3 1 52694E-03 



3.100 

-•31636442E-03 



3.200 

-• 1 2 515112E-03 



3.300 

•25396192E-03 



3.400 

• 1 1 672 168 E-03 



3.500 

-.20549376E-03 



3.600 

-•10984179E-03 



3.700 

• 1655 19 16E-03 



3.800 

• 10457756E-03 



3.900 

-•12943515E-03 



4.000 

-•98196927E-04 



4.100 

• 959 1088 1 E-04 



4.200 

•88 133242E-04 



4.300 

-.65782975E-04 



4.400 

-•73275802E-04 



4.500 

•40405606E-04 



4.600 

•53352386E-04 



4.700 

-.21199889E-04 



4.800 

-•28416106E-04 



4.900 

• 10 545189E-04 



5.000 

•13537622E-08 
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Indefinite integration and smoothing* Desired output and actual 
output. 


t 

Desired output 

,000 

-.3536779 3 E-00 

.100 

-.15340979E-00 

.200 

• 1 1219169E-00 

.300 

.37465 896 E-00 

.400 

,56 044806 E-00 

.500 

. 61228028E-00 

.600 

• 50686785E-00 

.700 

.26352543E-00 

.800 

-•59272709E-01 

.900 

-.37880758E-00 

1.000 

-.6 104846 5 E-00 

1.100 

-.69133740E-00 

1.200 

-.59788727E-00 

1.300 

-.35313085E-00 

1.400 

-•20540975E-01 

1.500 

• 3 1 327197E-00 

1.600 

•5623Q139E-00 

1.700 

.66451015E-00 

1.800 

. 59819582E-00 

1.900 

• 38682170E-00 

2.000 

•91 174926E-01 

2.100 

-.20860285E-00 

2.200 

-.43570072E-00 

2.300 

5374080 IE-00 

2.400 

-.49 8 10027 E-00 

2.500 

-.34105079E-00 

2.600 

-.11908317E-00 

2.700 

• 10262141E-00 

2.800 

•26607147E-00 

2.900 

.33641861E-00 

3.000 

.30975637E-00 

3.100 

.21090859E-00 

3.200 

.82 557013E-01 

3.300 

-.30313693E-01 

3.400 

-.95377830E-01 

3.500 

-.10249805E-00 

3.600 

-.64953839E-01 

3.700 

-•12357850E-01 

3.800 

♦2 1974700E-01 

3.900 

• 14623757E-0 1 


t Actual output 


.000 

-.35144868E-00 

.100 

- • 14 844095 E-00 

.200 

.11 825555E-00 

.300 

• 38 025760E-00 

.400 

, 56484080E-00 

.500 

.61 527389E-00 

.600 

, 50 8 14904E-00 

.700 

.26264061E-00 

.800 

-.62261615E-01 

.900 

-.38298996E-00 

1.000 

-.61 469266 E-00 

1.100 

- .6949 88 35E -00 

1.200 

-.60102638E-00 

1.300 

-.35568815E-00 

1.400 

-.2 1 8 68438 E -01 

1.500 

.31 391495E-00 

1.600 

, 56495094E-00 

1.700 

• 66844175E-00 

1.800 

• 60 263503E-00 

1.900 

• 39 141 123E-00 

2.000 

• 9 5 579163E-0 1 

2.100 

-.20530551E-00 

2.200 

-.43481 667E -00 

2.300 

-. 5 3 966433 E-00 

2.400 

-.50302054E-00 

2.500 

-.34737171E-00 

2.600 

-.12560173E-00 

2.700 

, 968 19607E-0 1 

2.800 

.26204345E-00 

2.900 

.33549022E— 00 

3.000 

• 3 1 278276E-00 

3.100 

.2 1739711E-00 

3.200 

.90734710E-01 

3.300 

-.22491099E-01 

3.400 

-.89373466E-01 

3.500 

-. 992343 14E -01 

3.600 

-•65210438E-01 

3.700 

-.16574201E-01 

3.800 

. 14440005E-01 

3.900 

• 57324858E-02 



BIBLIOGRAPHY 


1 . Titchmhrsh, Edward Charles, Introduction to the Theory of Fourier 

Integrals , The Clarendon Press, Oxford, 1937 . 

2. Tolstov, Georgi P. , Fourier Series , Prentice-Hall, Englewood Cliffs, 

N. J., 1962. 

3* Lighthill, M. , Introduction to Fourier Analysis and Generalized 
Functions , Cambridge [Eng.], University Press, 1958. 

4. Friedman, Avner, Generalized Functions and Partial Differential 

Equations , Prentice -Hall, Englewood Cliffs, N. J., 1963 . 

5 . Papoulis, Athanasios, The Fourier Integral and Its Applications , 

McGraw-Hill, New York, 1962. 

6 . Martin, Marcel A. , "Frequency Domain Applications in Data 

Processing", GE Technical Information Series, 57SD340, May 1957- 

7* Martin, Marcel A. , "Digital Filters for Data Processing, GE 
Technical Information Series," 62SD484, October 1962. 

8. Wiley, C. R. , Jr., Advanced Engineering Mathematics , McGraw-Hill, 

New York, i 960 . 

9* Morse, Phillip M. and Feshbach, Herman, Methods of Theoretical 
Physics , McGraw-Hill, New York, 1953* 

10 . Carslaw, Horatio Scott, Introduction to the Theory of Fourier ' s 

Series and Integrals , Dover, New York, 1930. 

11. Ormsby, Joseph F. A., "Design of Numerical Filters with Applications 

to Missile Data Processing," Journal of the Association for 
Computing Machinery, Vol. 8, No. 3 > July 1961. 

12. Graham, Ronald J. , "Determination and Analysis of Numerical Smoothing 

Weights," NASA Technical Report, R-179, 1963 . 


121 



"Remaiqur Concernint 1' integration 


13. Stiefal, E. and Ritishauser, H. , 

Numerique," Comptes Rendus, Academie Des Sciences 252, 
March-April 1961, pp. 1899-1900. 

14. Anders, E. B. , "Design of Numerical filters for more than one 

variable," unpublished. 

15. Friant, R. J., Jr., "Practical Digital Data Smoothing," GE 

Technical Information Series, R58EMH29* June 1958* 

16. Fleck, J. T. and Fryer, W. D. , "An Exploration of Numerical 

Filtering Techniques," Cornell Aeronautical Laboratory, Inc., 
Report No. XA-869-P-1* 


122 


* 


APPENDIX A 
CONSTRAINTS 


In order to develop constraints on the weights such that 
the recovered transfer function H has an exact fit at some specified 
frequency r we need to consider two separate cases. The first is 


N 


f * 


when H is of the form H(r) = h + 2 ) h cos 2nnr, r = — = 

n / n 7 


L -1 n 

n=l 


f 2*f * 

s s 


N 

The second is when H(r) is of the form H(r) = 2i^h n sin 2«nr. 

n=l 


A.l Constraints at one point 


N 


Case I. Suppose H(r) = h + 2 > h cos 2«nr. 
o n 


n=l 


then 


N 

H' (r) = nh^sin 2nnr. 

n=l 

We wish to impose the following constraints: 
H(r) = F(r) , 


X • 6 • y 


H'(r) = F'(r), 


N 

h Q + 2 h^cos 2itnr - F(r) = 0, 
n=l 
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4* ) nH sin 2«nr + F'(r) = 0. 

L “ 

n=l 


In order to minimize the error between H and H under the above con- 
straints we define 


R = / [H(r) - H(r)] 2 dr + a[4rt \ Hh^sin 2*nF + F'(r)]. 


Since 


h = F(r) - 2 ) h cos 2«nr, 

o / i n 

1 n=1 

2 N N 

R = / [F(?) + 2 ) h n (cos 2nnr - cos 2nnr) - h Q - 2 J h n cos 2«nr] dr 

0 n=l n=l 


+ a [bit \ nE n sin 2jtriF + F'(y)]. 


- ? ■■■ s 2 J [F(r) + 2\ h n (cos 2*nr - cos 2rtnr) - h Q - 2^h n cos 2itnr] 
d \ 0 n=l n=l 


[cos 2*kr - cos 2*kr]dr + a[4«ksin 2«kr]. 


Let = 0, k = 1, . . .,N. 

dh k 


h. A h 

i F(t)cos 2rtkT + "§ + ) h^cos 2nnT cos 2*kT + -^cos 2rtkT 


a[nksin 2itkr] 
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N 


^(l^-h^) + [ ^h n cos 2rtn? - -- ^- -]cos 2«kr 
n=l 

h 


+ -gCos 2rtkr = - a(rtksin 2«kr]. 


Let 8 = (h -h ), then 
o o 


(l^-h^) = 8 cos 2itkr - 2a[nksin 2nkr]. 

Multiply (A.l) hy (2 cos 2rtkr). Summing from 1 to N gives 

N N N 

2 


(A.l) 


X<W 


2 ) (h,_-h,Jcos 2rtkY = 2 ^ cos 2itkr " to / rtkcos 2«kr sin 2«kr, 
k=l k=i k=i 


adding (h o -h Q ) to both sides gives 


W 

<VV + 2 MW cos 2nkr 
k=l 

N N 

= ^ + 2 S ) cos^ 2«kr - to^rtkcos 2nkr sin 2«kr. 

k=l k=l 


Let 


N 


A^ = h Q + 2 y h n cos 2nnr - F(r) . 

n=i 


Hence 


N 


A, = (h -h ) + 2 ) (h -h )cos 2«nr, 
ML o o / , n n 7 

n=l 


so 


N 


N 


A. = 4a ) jtncos 2jtnr sin 2nn? - $ - 2 cos 2«nr. 
n=l n=l 


125 



Now multiply (A.l) by 2ksin 2«kr. Summing from 1 to N gives 

N N N 

2 ^ k ^ h k"\^ sin 2,tk * = 2 $/^ kcos 2 itkr sin 2 itkr - sk 2 sin 2 2 rtkr. 

3t=l k=l k=i 


Hence 


1M 

= -4^nh n sin 2 rtnr - F'(r). 


IM 

Ag = 4it ^n(h n -h n )sin 2 *nr. 


Ag = 4jt S^kcos 2 rtkr sin 2 rtkr" - 8 it 2 a ^ k 2 sin 2 27 tkr. 


Let = 2 y cos 2rtkr, 


q 2 = 4„v kcos 2 «kr sin 2 nkr, 


oP 2 2 _ 

Q 3 ■ 8 « \ k sin 2 rtkr. 


\ - QgG - (1+^) > 


A 2 = -e 3 a + 02 fi 
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Solving we find that 


o 2, ^ — > 

Q^-(1+Q 1 )Q 3 

- 

4 -^% 

Therefore the constrained weights are 

h = h + £ , 
o o 


= h^ + $ cos 2*kr - a2«ksin 2nkr, k > 1, 

where £ and a are as defined in (A. 2) and (A. 3). 

N 

Case II . Suppose H(r) = 2i^h n sin 2«nr, 

n=l 


then 


N 

H'(r) = 4ni 

n=l 

We wish to impose the following constraints 


Z nh cos 2jtnr. 

n 


i.e. 


> 


and 


H(r) ? F(r), 

H'(r) = F'(r), 

N 

2^h n sin 2*nr - = 0, 

n=l 

N 

bn ) nh cos 2itnx - ■ - r = 0 . 

/_, n 1 

n=l 


(A. 2 ) 


(A. 3) 
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In order to minimize the error between H and H under the above 

conditions we define 

1 
2 


R = j^[H(r) - H(r)j dr + a[4it^nh n cos 2nnr - — ]. 

0 n=l 


N 


Since 


F( 


h-. 


n=2 


sin 2flnr 


2sin 2icr 


1 

2 


2 N N 

R = f [ [ F i r ^ - 2 V h sin 2itnr ] + 2 \ h sin 2«nr 

w e4 n O-J. 1 L - I n l—l n 

0 sin dnT n=2 n=2 


™ _ lill 


2 N 

- 2 ) h sin 2«nr] dr + a[4n ) nh cos 2«nr - - - ] 

/_, n L, n 1 

n=l n=l 

1 

2 N N 

f [ sin_2«£ [ . 2 y h sin 2*nr] + 2 V i^sin 2*nr 

d \ 0 sin 2 ” r n'h ” n=2 " 


N 

-2 ) h sin 2«nr][ ~ 2sin — s ^ n 2Ttkr + 2sin 2*kr]dr + 4jtOfecos 2«kr. 

L—l n sin 2«r 

n=l 


Let = 0, k = 2, . 

K 

1 

2 

Now f sin^lM- 
1 sin 2rtr 1 


. .,N. 


N 

_Vr j ~ -u . _ . sin 2«r sin 2itkr-i , 

2 ) h sin 2«nrJ[sin 2nkr - — Jdr 

£-± n sin 2irr 

n=2 

N 

= . £©■ Bln! >*kT + sin. a&r y jj sin 2l(n? 

1 sin 2 2«r 2sin 2«r (-L n 
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1 
2 N 


Also / 1®. -h^lsin 2jtnr[sin 2«kr - 


sin 2«r sin 2«kr-, l r r- 
- ] dr = 


0 n=2 

1 
2 


sin 2nr 


r [h sin^ 2nr] s ^ - dr = 4 ^^5 . 

J 1 sin 2nr 4 


sin 2jtr 


Hence 


» _ h 

i [2 V h sin 2*nr - IM] + ± ft-h ] ♦ 1 

4 n 4^ sin 2*r 4 K k 4 sin 2*r 

n=<£ 


a 


- 


cos 2rtkr 


( h l~hi) — 2]X ^ k + (h^-h^) = - a«kcos 2jtkr. 


sin 2nr 


Let 8 = then 


h -h = $ — 2n ^ - ankcos 2nkr. (A. 4) 

k sin 2*r 

Multiplying (A. 4) by 2sin 2rtkr and summing from 2 to N gives 

N 1 2 — 

2-V (h, -h )sin 2«kr = $ 2 \ 1 - 02 \ k cos 2nkr sin 2«kr. 

£*> k^2 sin 2 * r * k=2 


Adding 2(h^-h^) sin 2itr to both sides yields 


_N_ N 

2 \ (h, -h ) sin 2jtkr = 2 §\ s ^ n ■ - nkr -2a \ kcos 2«kr sin 2nkr. 
-U * k L s in 2*F .U 


k=l 


k=l sin 2flr k=2 


N 


Let = 2^h k sin 2«kr - F(r). 
k=l 
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Since 


or 


N 

F(r) = 2 ^ h^sin 2nkr, 
k=l 

N 

A = 2 I ( V\> sin 2nkr, 
k=l 


N 


N 


2 _ 
sin 2«kr 


A, = 20t ) kcos 2«kr sin 2nkr - 2 $ _ 

- 1 - l-l ,^-i sin 2«r 

k=2 k=l 

Multiplying (A. 4) by 4nkcos 2*kr and summing from 2 to N gives 

N 


4rt 


I'VV 


kcos 2«kr = 8 4« 


N N 

ksin 2rtkr cos 2wkr 


sin 2nr 


-hoc ) « 


2 , 2 2 

k cos 


k=l k=l k=2 

adding 4«(h^-h^)cos 2«r to both sides of the above equation gives 


N 


4« (h^-h^kcos 2nkr = 8 4* ^ 
k=l k=l 


N _ N 

kcos 2«kr sin 2«kr ^ V 1 rt ^.^ cos ^ 


sin 2«r 


L, 

k=2 


N 

Let Ag = 4*^101^008 2«kr - F'(r). 
k=l 

N 

Since F'(r) = 4n^kh k cos 2*kr. 
k=l 


N 

a 2 = ' t ”I ( V\ ) kcos 2nkr. 
k=l 


Hence 


N _ N 

g Y kcos 2nkr sin 2* kr , 2 V ,_2 2 

sin 2itr 


a ), C \ tllUVX SJLLI c:hjxx c. \ . c. c. , — 

A^ = 4jt 6 ) 4an > k cos 2*kr. 


l 1 

k=2 


2«kr, 


2nkr. 
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Let 


N 


Ql = 2 


k=l 


2 — 

sin 2rtkr 

sin 2jtr 


N 

Q 2 = 2^kcos 2nkr sin 2nkr, 
k=l 


S’ 


sin 2jtr 


Then 


and 


N 

Q,^ = 4* 2 ^ k 2 cos 2 2«kr. 
k=l 


A^ = - cos 2itr sin 2itr) - $ Q^, 


Ag = -a(Q,^-4rt 2 cos 2 2«r) + $ Q_ 
Solving for 8 and a we find that 


2 2 

cos 2nY) + A^Wg-cos sin 2«r) 


8 = 


A p 

Q^(Q 2 -cos 2*r sin 2jtr) - OgCO^-ta cos 2«r) 


a = 


Aj Q ^ " ^2 


Q^^-cos 2ltr sin 2jtr ) -Q 2 (Qn."^ rt2cos2 2«r) ) 
Therefore the constrained weights are 


h l = h l + $ 


vT - v> . e* sin 2«rk „ . „ , — . _ _ 

n^ - ^ + S — - ank cos 2rtkr, k > 2. 

sin 2jtr 


(A.5) 


(A .6) 
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APPENDIX B 

ROMBERG'S METHOD OF NUMERICAL INTEGRATION 

b 

Suppose we wish to calculate J = J g(x)dx. First we define 


2 k -l 


T q k = 2~ k [ | g(a) g(a+ n) + | g(b)], 

n=i 2 

which are the results of applying the trapezoidal rule with a 
partitioning of the interval [a,b] into 2 parts. With these values 


we form 


JL . 

m,k 


k^T -T 

* m-1, k+1 m-l,k 


U m -1 


m=l,2, k=0, 1,2, and arrange them in the triangular table 

T 

0,0 

T T 

0,1 1,0 

T T T 

o ,2 1,1 2,0 


n T T T 

o,3 1,2 X 2,l A 3 ,o 


Then J = (b-a) T . See [133 • 

Ulj o 

k 

If the interval [a,b] is divided into 2 subintervals, each of 

2k+2 

length h, it has been shown that J - (b-ajT^ Q = 0(h ). 

We note that, assuming that a sufficiently large number of 
functional values are available in the interval, a very close 
approximation of the integral can be obtained. 
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